

Sign up
Here's how it works:
  1. Anybody can ask a question
  2. Anybody can answer

The formal problem

The question I am about to ask is applicable to many scientific and engineering problems. In short, I am bout to ask for pointers and practical advice for Mathematica beginners who would like to use Mathematica to find a solution to a multi-dimensional non-linear constraint optimization problem, i.e. a problem of the form:

Find $\vec{\mathbf{x}}$ so that $\mathbf{F}\left(\vec{\mathbf{x}}\right) = \vec{\mathbf{0}} $. To expand this expression:

Find a value of $\vec{\mathbf{x}} = \left| \begin{array}{c} x_{1} \\ x_{2} \\ \vdots \\ x_{n} \end{array} \right|$ so that $\begin{cases} f_1\left(x_0, x_1, \ldots, x_n\right) = 0 \\ f_2\left(x_0, x_1, \ldots, x_n\right) = 0 \\ \mbox{ }\vdots& \\ f_m\left(x_0, x_1, \ldots, x_n\right) = 0 \end{cases}$

Typically, we are looking for $\vec{\mathbf{x}} \in \Re^n$ (i.e. ignore complex solutions in most engineering and scientific disciplines), we have $n \neq m$, and $n$ and $m$ are of single-digit to two-digit numbers. The functions $f_i$ are arbitrary, non-linear and only defined on subspaces of $\Re^n$. There may be additional inequality constraints.

The practical problem

I would like to know all the tools, functions and tricks experienced Mathematica users like to use to approach this kind of problems. I am not asking for a one-command magic "solve that for me" function, as I know that function optimization is in general very hard, and a combination of formal analysis, plotting and numerical exploration are what guides one to a solution, if a solution exists. What I am asking is how Mathematica can help this exploration.


Let's take a real-world example in Mathematica syntax, taken form a problem on which I am currently working. I believe this particular set of equations has no real solution (so don't try to actually work on it!) and could be better formulated, but that makes it an even better example.

The system looks like that (you don't need to download the file, it's not very interesting):

eq1 = Dhtmp==Chbar+Ghtmp
eq2 = Dfstartmp==Cfstarbar+Gfstartmp
(* ... *)
eq31 = Lmbar==((Dmtmp+Dmstartmp*((1-n)/n))/Wbar)*(1-zprimebar)
eq32 = Lmstarbar==((Dmtmp+Dmstartmp*((1-n)/n))/Wstarbar)*zprimebar*(n/(1-n))
eq33 = Wstarbar==((t/tstar)*(((1-zprimebar)/zprimebar)^(1/theta)))*Wbar
eq34 = zeta==(Ptmp*(Ctmp^phi))/(Pstartmp*(Cstartmp^phi))
eq35 = Wbar==Phbar/(Pmbar*ahbar)
eq36 = Wstarbar==Pfstarbar/(Pmstarbar*afstarbar)
eq37 = ibar==(1/beta)-1

This is the raw set of equations, in that case, a copy-paste from a Matlab file followed by a bit of search/replace to remove underscores, change implicit functions to explicit, label them, etc... With that, some pre-determined values or intervals are given for some of the variables, like ibar = 1, etc...

As a beginner Mathematica user, the first "trick" that I've learnt that lets me look at the system (thanks to MMA.SE user Simon) is to define a set of pretty-printing "replacements" for the ugly matlab variable names:

MakeBoxes[afstarbar, form_] := 
   InterpretationBox[OverscriptBox[SuperscriptBox["af", "*"], "_"], afstarbar]
MakeBoxes[alpha, form_] := InterpretationBox["\[Alpha]", alpha]
(* ... there are 53 of those; I paste them before the equations. *)

In Matlab, typically, I try to plot each equation individually, holding all but one or two variables constant. It helps me get a better idea of the search intervals of the free variables. Then I feed simplified versions of the system to various solvers (fsolve for simple stuff, or I write some simulated annealing-like one-time-use solver, etc...), and tune the system until I get a solution.

Now, I want to give Mathematica a try, but as a beginner, I know that it's going to be a lot of fumbling around, because I did not find a good list of pointers on what Mathematica users use in real-world to help themselves solve non-linear multi-dimensional constraint optimization problems.

I think that the Mathematica Stackexchange website is where such a guide should be found.

share|improve this question
MATLAB's fsolve() and Mathematica's FindRoot[] are essentially equivalent, I've seen. NMinimize[] implements simulated annealing and a bunch of other global optimization methods; you might want to look into them too. – J. M. Jan 26 '12 at 4:45
Thanks @j-m, I'll look more closely at NMinimize and FindRoot. I also read bout NSolve, I'll look at that one too when I have access to Mathematica. The best would be to find an online tutorial where one shows how he used Mathematica to find the solution to a non-trivial example. You know, step-by-step, with the errors and trials. – agravier Jan 26 '12 at 5:02
NSolve[] is a bit specialized for algebraic equations; if you have transcendentals involved, you'd want to use FindRoot[] instead. – J. M. Jan 26 '12 at 5:04

Your system is very large and what you're asking is quite general. So I guess only collection of general tips or links can help for now. If you provide a more specific smaller example - it would be easier to address. Generally there is an extensive documentation with many examples on "Constrained Optimization":

An example of a flexible approachs to a solution of a nonlinear system was already posted here:

Updating Wagon's FindAllCrossings2D[] function

List of some relevant interactive Demonstrations with free source code that solve a nonlinear system:

share|improve this answer
Thanks for your answer. Yes, it's a general question, I'm looking for general advice and worked out real examples. I already noticed and I like a lot the function that uses the graph to guide to a solution. TBH, I don't really need to solve such large system right now, but I think that it would be beneficial to all noobs to collect resources and links in one thread, that's why I asked such a general question. – agravier Jan 28 '12 at 4:19

Your Answer


