Using the Solve code in your most recent message:
Each solution is enclosed in {} and all the solutions are a list enclosed in another layer of {} so I can use
Length[sol]
to find the length of that list and it tells me there are 91 different solutions.
Next, there is nothing in the Mathematica code of your most recent message that says a2x must be greater than or equal to zero. If you want to include that condition then this
sol=Solve[{-2a1x+10/(a1x+a1y-10)+g1==0,
-2a1y+10/(a1x+a1y-10)+g2==0,
9/10(-2a2x(1-ix/40)+10/(a2x+a2y-10))+g3==0,
9/10(-2a2y(1-ix/40)+10/(a2x+a2y-10))+g4==0,
-1/2+9/10(a2x^2/40+a2y^2/40)+g5==0,-1/2+g6==0,
g1*-a1x==0,g2*-a1y==0,g3*-a2x==0,g4*-a2y==0,
g5*-ix==0,g6*-iy==0,g7*(ix-24)==0,g8*(iy-24)==0,
a2x>=0},
{a1x,a1y,a2x,a2y,ix,iy,g1,g2,g3,g4,g5,g6,g7,g8}];
Length[sol]
reports that there are now only 63 solutions.
If I go back to an earlier message of yours and use all the conditions that it did and hope that this is the right thing to do then
sol=Solve[{-2a1x+10/(a1x+a1y-10)+g1==0,
-2a1y+10/(a1x+a1y-10)+g2==0,
9/10(-2a2x(1-ix/40)+10/(a2x+a2y-10))+g3==0,
9/10(-2a2y(1-ix/40)+10/(a2x+a2y-10))+g4==0,
-1/2+9/10(a2x^2/40+a2y^2/40)+g5==0,-1/2+g6==0,
g1*-a1x==0,g2*-a1y==0,g3*-a2x==0,g4*-a2y==0,
g5*-ix==0,g6*-iy==0,g7*(ix-24)==0,g8*(iy-24)==0,
a1x>=0,a1y>=0,a2x>=0,a2y>=0,0<=ix<=24,0<=iy<=24,
g1>=0,g2>=0,g3>=0,g4>=0,g5>=0,g6>=0,g7>=0,g8>=0},
{a1x,a1y,a2x,a2y,ix,iy,g1,g2,g3,g4,g5,g6,g7,g8}];
Length[sol]
tells me that there are only two solutions and thus
sol
shows me
{{a1x->(5 + Sqrt[35])/2, a1y->(5 + 5*(5 + Sqrt[35]))/(5 + Sqrt[35]),
a2x->0, a2y->0, ix->0, iy->0, g1->0, g2->0, g3->9/10, g4->9/10, g5->1/2, g6->1/2, g7->0, g8->0},
{a1x->0, a1y->0,
a2x->0, a2y->0, ix->0, iy->0, g1->1, g2->1, g3->9/10, g4->9/10, g5->1/2, g6->1/2, g7->0, g8->0}}
as the only two solutions. Note: It might not be important, but if you Simplify that value for a1y in the first solution it turns out to be the same as the value for a1x, it just looks different because some functions don't take the extra time to do simplification that you may not need or want done.
If on the other hand the only constraint you want to put on the solutions is that ix==24 then include that, as the only additional equation.
sol=Solve[{-2a1x+10/(a1x+a1y-10)+g1==0,
-2a1y+10/(a1x+a1y-10)+g2==0,
9/10(-2a2x(1-ix/40)+10/(a2x+a2y-10))+g3==0,
9/10(-2a2y(1-ix/40)+10/(a2x+a2y-10))+g4==0,
-1/2+9/10(a2x^2/40+a2y^2/40)+g5==0,-1/2+g6==0,
g1*-a1x==0,g2*-a1y==0,g3*-a2x==0,g4*-a2y==0,
g5*-ix==0,g6*-iy==0,g7*(ix-24)==0,g8*(iy-24)==0,
ix==24},
{a1x,a1y,a2x,a2y,ix,iy,g1,g2,g3,g4,g5,g6,g7,g8}];
Length[sol]
and that returns zero, because sol is {} because apparently there are no solutions when ix==24.
So, is what I've shown you thus far the solution to your problem? Or at least showing you the method that I use to attack problems like this using the Mathematica way of thinking?
As for "simultaneously optimize two objective functions that are dependent on each other", I don't have any confidence that I know how to translate that phrase into a line or two of Mathematica code. If you can show those two functions written in Mathematica notation and you can describe what your optimization model is then perhaps I can try to understand that and see if I can translate it into a few lines of Mathematica code.
Is what I am doing with this helping you at all?
Is there any way you could do what you are trying to do entirely within Matlab? Is there any chance that working within one tool and language and way of thinking might be easier for you than trying to translate back and forth between two tools and languages and ways of thinking?