---

11.115 Sample ‘solve’ ¶

Example of solve for root finding.

MGL code:

zrange 0 1
new x 20 30 '(x+2)/3*cos(pi*y)'
new y 20 30 '(x+2)/3*sin(pi*y)'
new z 20 30 'exp(-6*x^2-2*sin(pi*y)^2)'

subplot 2 1 0:title 'Cartesian space':rotate 30 -40
axis 'xyzU':box
xlabel 'x':ylabel 'y'
origin 1 1:grid 'xy'
mesh x y z

# section along 'x' direction
solve u x 0.5 'x'
var v u.nx 0 1
evaluate yy y u v
evaluate xx x u v
evaluate zz z u v
plot xx yy zz 'k2o'

# 1st section along 'y' direction
solve u1 x -0.5 'y'
var v1 u1.nx 0 1
evaluate yy y v1 u1
evaluate xx x v1 u1
evaluate zz z v1 u1
plot xx yy zz 'b2^'

# 2nd section along 'y' direction
solve u2 x -0.5 'y' u1
evaluate yy y v1 u2
evaluate xx x v1 u2
evaluate zz z v1 u2
plot xx yy zz 'r2v'

subplot 2 1 1:title 'Accompanied space'
ranges 0 1 0 1:origin 0 0
axis:box:xlabel 'i':ylabel 'j':grid2 z 'h'

plot u v 'k2o':line 0.4 0.5 0.8 0.5 'kA'
plot v1 u1 'b2^':line 0.5 0.15 0.5 0.3 'bA'
plot v1 u2 'r2v':line 0.5 0.7 0.5 0.85 'rA'

C++ code:

void smgl_solve(mglGraph *gr)	// solve and evaluate
{
	gr->SetRange('z',0,1);
	mglData x(20,30), y(20,30), z(20,30), xx,yy,zz;
	gr->Fill(x,"(x+2)/3*cos(pi*y)");
	gr->Fill(y,"(x+2)/3*sin(pi*y)");
	gr->Fill(z,"exp(-6*x^2-2*sin(pi*y)^2)");

	gr->SubPlot(2,1,0);	gr->Title("Cartesian space");	gr->Rotate(30,-40);
	gr->Axis("xyzU");	gr->Box();	gr->Label('x',"x");	gr->Label('y',"y");
	gr->SetOrigin(1,1);	gr->Grid("xy");
	gr->Mesh(x,y,z);

	// section along 'x' direction
	mglData u = x.Solve(0.5,'x');
	mglData v(u.nx);	v.Fill(0,1);
	xx = x.Evaluate(u,v);	yy = y.Evaluate(u,v);	zz = z.Evaluate(u,v);
	gr->Plot(xx,yy,zz,"k2o");

	// 1st section along 'y' direction
	mglData u1 = x.Solve(-0.5,'y');
	mglData v1(u1.nx);	v1.Fill(0,1);
	xx = x.Evaluate(v1,u1);	yy = y.Evaluate(v1,u1);	zz = z.Evaluate(v1,u1);
	gr->Plot(xx,yy,zz,"b2^");

	// 2nd section along 'y' direction
	mglData u2 = x.Solve(-0.5,'y',u1);
	xx = x.Evaluate(v1,u2);	yy = y.Evaluate(v1,u2);	zz = z.Evaluate(v1,u2);
	gr->Plot(xx,yy,zz,"r2v");

	gr->SubPlot(2,1,1);	gr->Title("Accompanied space");
	gr->SetRanges(0,1,0,1);	gr->SetOrigin(0,0);
	gr->Axis();	gr->Box();	gr->Label('x',"i");	gr->Label('y',"j");
	gr->Grid(z,"h");

	gr->Plot(u,v,"k2o");	gr->Line(mglPoint(0.4,0.5),mglPoint(0.8,0.5),"kA");
	gr->Plot(v1,u1,"b2^");	gr->Line(mglPoint(0.5,0.15),mglPoint(0.5,0.3),"bA");
	gr->Plot(v1,u2,"r2v");	gr->Line(mglPoint(0.5,0.7),mglPoint(0.5,0.85),"rA");
}