A pragmatic approach is just to use off-the-shelf black-box mathematical optimization algorithms.
In particular, define parameters $\theta_1,\theta_2,\theta_3$, representing the angle from the horizontal of $D_1,D_2,D_3$, respectively. Here each $\theta_i$ is in the range $[0,2\pi)$. The locations of $D_1,D_2,D_3$ are completely determined by the value of $\theta_1,\theta_2,\theta_3$. Therefore, we can express the value of $F$ as a function of $\theta_1,\theta_2,\theta_3$: let's write this as $F^*(\theta_1,\theta_2,\theta_3)$. It's tedious but not too hard to write down an explicit expression for $F^*$.
The problem becomes: given an expression for $F^*$, find $\theta_1,\theta_2,\theta_3$ that maximize $F^*(\theta_1,\theta_2,\theta_3)$, subject to $0 \le \theta_1,\theta_2,\theta_3 < 2\pi$. This is a mathematical optimization problem on a three-dimensional objective function $F^* : [0,2\pi)^3 \to \mathbb{R}_{\ge 0}$. There are standard algorithms for maximizing such a function, e.g., using gradient descent or any number of other methods. You can try applying any of them.
In practice, I suggest you run gradient descent, starting from an initial point (an initial assignment to $\theta_1,\theta_2,\theta_3$) chosen at random. Run gradient descent until convergence. This will give you a local maximum. Repeat this 1000 times, with 1000 different randomly chosen starting values, and take the best solution found in any of these iterations. While there are no guarantees, I expect this to work well. Generally speaking, these methods tend to work very well when the number of dimensions is small, and in this case, you have only 3 variables and 3 dimensions, so I expect they'll work well. In particular, I expect the function $F^*$ to be piecewise smooth, with only a small finite number of pieces, so this procedure will probably ensure that you have at least one starting point in each piece and thus that gradient descent finds the maximum within each piece, ensuring that you find the global maximum with high probability.
If you wanted to do a more careful analysis, you could probably work out analytically the boundaries of each piece (i.e., decompose the three-dimensional space $[0,2\pi)^3$ into a finite union of regions within which $F^*$ is continuous and differentiable) and then run one iteration of gradient descent per piece. However, that sounds tedious, and it's probably easier to just slam the thing into an off-the-shelf optimization tool and let the computer do the grunt work for you. Doing the delicate analysis might be important if you want a procedure that is guaranteed correct in all cases (e.g., long and skinny triangles, etc.), but if you just want something that'll probably work without working too hard, I suggest you spare yourself the headache and let the computer do the crunching.
