--- a/example/Elasticity/HyperElasticity.m	Mon Aug 12 22:13:37 2013 +0200
+++ b/example/Elasticity/HyperElasticity.m	Mon Aug 12 23:57:03 2013 +0200
@@ -10,26 +10,26 @@
 msho = msh3m_structured_mesh (x, y, z, 1, 1:6);
 mshd = Mesh (msho);
 
-global V  = HyperElasticity_FunctionSpace (mshd);
+V  = HyperElasticity_FunctionSpace (mshd);
 
-global bc2 = DirichletBC (V, @(x,y,z) [0; ...
+bc1 = DirichletBC (V, @(x,y,z) [0; 0; 0], 1);
+
+bc2 = DirichletBC (V, @(x,y,z) [0; ...
  0.5*(0.5 + (y - 0.5)*cos(pi/3) - (z-0.5)*sin(pi/3) - y);...
  0.5*(0.5 + (y - 0.5)*sin(pi/3) + (z-0.5)*cos(pi/3) - z)], 2);
 
-global bc1 = DirichletBC (V, @(x,y,z) [0; 0; 0], 1);
-
-global B = Expression ('B', @(x,y,z) [0.0; -0.5; 0.0]);
-global T = Expression ('T', @(x,y,z) [0.1; 0.0; 0.0]);
-global E = 10.0;
-global nu = 0.3;
-global mu = Expression ('mu', @(x,y,z) E./(2*(1+nu)));
-global lmbda = Expression ('lmbda', @(x,y,z) E*nu./((1+nu)*(1-2*nu)));
+B = Expression ('B', @(x,y,z) [0.0; -0.5; 0.0]);
+T = Expression ('T', @(x,y,z) [0.1; 0.0; 0.0]);
+E = 10.0;
+nu = 0.3;
+mu = Expression ('mu', @(x,y,z) E./(2*(1+nu)));
+lmbda = Expression ('lmbda', @(x,y,z) E*nu./((1+nu)*(1-2*nu)));
 u = Expression ('u', @(x,y,z) [0; 0; 0]);
 
 L = HyperElasticity_LinearForm (V, mu, lmbda, B, T, u);
 u0 = assemble (L, bc1, bc2);
-
-[x, fval, info] = fsolve (@f, u0, optimset ("jacobian", "on"));
+fs = @(xx) f (xx, V, bc1, bc2, B, T, mu, lmbda)
+[x, fval, info] = fsolve (fs, u0, optimset ("jacobian", "on"));
 
 func = Function ('u', V, x);
 fem_plot (func);

--- a/example/Elasticity/f.m	Mon Aug 12 22:13:37 2013 +0200
+++ b/example/Elasticity/f.m	Mon Aug 12 23:57:03 2013 +0200
@@ -1,12 +1,4 @@
-function [y, jac] = f (xx)
-
-  global V;
-  global mu;
-  global lmbda;
-  global T;
-  global B;
-  global bc1;
-  global bc2;
+function [y, jac] = f (xx, V, bc1, bc2, B, T, mu, lmbda)
 
   u = Function ('u', V, xx);