@Thomas Richard 
Your hint (?fnormal) is in line with a solution I got from another user: i.e. limiting the digits of the numeric evaluation of an integral.

Limiting digits most likely removed such terms.

I can see now fsolve having to deal with information on directions pointing away from the real valued solution...

Thanks for the quick clarification!

@Christian Wolinski This completes my picture. 

All your expressions are better than evaluation of the integral within fsolve (which I still think should work).

Your first expression is still the best in terms of a) computation efficiency and b) no restrictions on the search interval. I am using a squared form which removes the s shape of the function in the search interval.

g := alpha -> fsolve(alpha = EllipticF(sqrt(sin(x__0))/sqrt(sin(x__0) + 1), sqrt(-sin(x__0)^2 + 1)*I/(sin(x__0) - 1))^2/(1 - sin(x__0)),x__0)

In the meantime I have found another simple looking expression with the help of Maples Functionadvisor. Pretty cool (the advisor).

Int(1/sqrt(sin(x0) - sin(x)), x = 0 .. x0) = sqrt(2)*InverseJacobiCN(cos(x0)/(sin(x0) + 1), sqrt(2)*sqrt(sin(x0) + 1)/2).

Thanks again! Without your help I would have to live with numerical compromises.

@Carl Love
It works fine now

I can even access the output component wise! If it was for a Maple functionality I would suggest an access via indices where each index corresponds to either an operand or an operator (i.e. 1..9 for the example above). 

However, for the moment I am fine with the solutions as they are. Thank you all!

@Carl Love. I tried it but the output still contains the unwanted parenthesis.

Apart from that I like the extension to more than 3 operators and the inclusion of other operator types.

If I could later on refer in an easy way to components of the output, this would be perfect.
I vote up since it is close to a solution that I consider worth implementing in Maple.

Multiequation.mw

@acer  Thank you for your very detailed answer!

I have attached an example how I would use your solution 1). The input is hidden.

In case someone from MapleSoft reads this: Please take the attachment as a suggestion to improve Maple in the field of

- professional looking documents

- integration (computing the inverse with fsolve, simplifying integrals to expressions containing Maple functions)

Download inverse_function_with_fsolve_with_EllipticF.mw

@Christian Wolinski Thank you! Where did you get the expressions from? I was acutally expecting Maple to return such expressions.
I am interested, since I work in the meatime with the square of the expression (see my second attachment above). I guess that simliar expressions can be obtained.

Update: The simple formulas work much better. No reduced prescsision is required.

I pasted you code in my version and it worked.

Which version are you using?

I finally had to implement both proposals from @mmcdara and @acer. Calculating expressions containing elliptic integrals and their inverse seems to challenge sophisticated high-level commands like fsolve or plot. I guess that the default convergence criteria are simply set too high in this instance.  

On the other hand, Maple was able to integrate expressions containing a potential division by zero without limiting the precision of the integration or the integration range by a small fraction. Thumbs up for Maple.

My problem stemmed from a rather simple looking BVP containing an ODE that is similar to the simple nonlinear pendulum.

Those who are interested can find details in the attachment how this problem can be solved analytically and where I stopped calculating the inverse with fsolve.

Recalc_of_Belendez_et_al.mw

@acer I was looking for such an option within fsolve. Thank you!

@mmcdara I read in a paper that such problems are "easily" solveable in other Math packages. However they did not give more than a hint: use of numerical integration.  They probably  did soemthing simliar to acers solution.

@Carl Love Thanks for clarifing.
Further bellow in the worksheet there is the statement "Solve 2nd order ODE for phi(phi)". Should this be u(phi) then?

@Anthrazit If you can reproduce it: Do you see a dependence when you change the dpi settings (e.g. changing the monitor)?

I am asking because I have rendering issues when I switch monitors. It could be Windows related and may depend on compatibility settings of the program.

It seems that this test case is very challenging since the analytical solution oscillates and loops back to zero. 

A force applied at the free end of the beam instead of a momentum is probably more forgiving. Loops are not possible, and the leaver arm decreases with the deflection. Unless for this new test case someone can come up with analytical or PDE solutions, verified FEM code (which I do not have access to) could be used to confirm this guess.  

For planar models: Modifying the default settings of the flexible beam compoent in the example by setting the elastic coordinate in z direction to zero will speed up the calculation (about a factor of 10 in the example with the segmented beam). This will half the number of generalized coordinates the solver has to deal with. 

Likewise, increasing the elastic coordinates will increase the number of generalized coordinates.  

@vv Works. Thank you!