"Vectorized" fsolve

Hi,

I have a system of about 16 or so nonlinear equations that I would like to solve with fsolve. These equations contain two parameters, beta and delta, that I would like to vary and find the solution to the equations for each combination of these parameters. There are somewhere around 170 combinations of beta and delta(possibly more, depending on whether I need better resolution or not).

To be more clear, I want to make a carpet plot(same idea as the picture here http://www.vikingaero.com/images/range_carpet_plot.png ) of the variables in the nonlinear equations versus beta and delta. In the picture, delta C sub D0 and delta W sub T0 would be beta and delta and R would be one of the variables that I'd be solving for in the non linear equations.

Is there any way that I can "vectorize" with fsolve, or do I have to loop through all my beta and delta combinations and call fsolve for each combination?

[Sean de Wolski]

How slow is the for-loop and how important is speed?

Keep it simple and stupid.

[Matt]

I don't know how slow it is for sure yet. I'm still debugging the solution and am having issues with convergence at some points. I was planning to go the simple and stupid way until I got it working for sure and then try to speed it up later. I figured that it would take awhile to get a response on here though(apparently I was wrong :) ), so I asked early.

Speed will be pretty important when all is said and done. I'd like to be able to do studies on how the many other parameters that go into these equations affect the results. To put it in perspective, this is a vehicle dynamics model and I would like to be able to run this on my computer at the track so that I can talk to the driver about how the car feels at various points on the track, test my proposed improvements to the car setup in this model and then decide on what changes to make to the car setup.

Right now, when it gets to a set of inputs for which it converges, it pumps out the "equation solved" messages faster than I can read them.

Answer 1

[Sargondjani]

first of all, if you loop through R and use the solution of the last loop as starting value, that will already help quite abit (less iterations needed compared to using same starting values)

but you could do without the loop by just specifiying alll 170 equations + 2 x 170 input variables and solve simultaniously. For speed you absolutely need to let matlab know that all equations are independent, so you have to either set the JacobPattern (for finite difference, otherwise it will take ages) or supply the analytical jacobian.

let us know if you have problems settings the JacobPattern (post some of your code)

[Matt]

Sorry, I should have been clearer. I have a set of 16 non-linear equations in 16 unknowns(the equations contain more than one unknown each) that I need to solve for each beta and delta combination. Even if I tried hard coding all 16 * 170 = 2720 equations, they would not be independent. Equations 1-16, 17-32, 33-48 and so on would be dependent on each other.

Also, the number of each beta and delta is not fixed. Depending on other input variables, I may find I need a wider or smaller range of beta or delta or a finer resolution.

Many elements of the equations are either lookup tables, fitted polynomials or other curve fits, so I can't supply an analytical jacobian.

I'll use a nearby solution value as a starting value for the next solution for now.

[Sargondjani]

i personally would not try to solve all 2720 equations simultaniously, as solving 16 equations for 16 unknows is already quite intensive... you could try multiple workers though (look up documentation on 'parfor')

however, you should definitely supply the "jacobpattern" for the finite differences for each 16 equations. this will be crucial to improve computation speed

[Matt]

I didn't even think of looking at parallel. That could be pretty good - I could give each worker a line of constant beta or a line of constant delta to work on. I'm running an i7-2600K and 2 GTX 560s in SLI, so I'll definitely give that a look to try and find more speed once I have the solution working 100%.

Doesn't fsolve automatically compute the Jacobian with finite differences? How would me supplying it be faster?

[Sargondjani]

yes fsolve does automatically compute the jacobian. however it assumes that it is a full matrix, which it will not be unless all 16 unknowns show up in every equation. and since you can know in advance which values in the jacobian will be zero, you should provide the zeros otherwise matlab will try to estimate them (use sparse matrix)....

for example IF all your equations would be independent you would just have to calculate the diagonal of your jacobian, instead of a full 16x16 matrix. see the point?