% Space dimensions
fem.sdim={'x','y','z'};

% Geometry
g1=block3(1,1,0.10000000000000001,'corner',zeros(1,3),[0 0 1],0);
g2=block3(1,1,0.40000000000000002,'corner',zeros(1,3),[0 0 1],0);
g3=block3(1,1,0.20000000000000001,'corner',zeros(1,3),[0 0 1],0);

clear s f c p
objs={g3};
name={'BLK1'};
s.objs=objs;
s.name=name;

objs={};
name={};
f.objs=objs;
f.name=name;

objs={};
name={};
c.objs=objs;
c.name=name;

objs={};
name={};
p.objs=objs;
p.name=name;

drawstruct=struct('s',s,'f',f,'c',c,'p',p);
fem.draw=drawstruct;
fem.geom=geomcsg(fem);

% Initialize mesh
fem.mesh=meshinit(fem,...
	'Out',    {'mesh'},...
	'jiggle', 'on',...
	'Hcurve', 0.40000000000000002,...
	'Hcutoff',0.01,...
	'Hgrad',  1.3999999999999999,...
	'Hmaxfact',1);

% Refine mesh
fem.mesh=meshrefine(fem,...
	'out',    {'mesh'},...
	'rmethod','longest');

% Refine mesh
fem.mesh=meshrefine(fem,...
	'out',    {'mesh'},...
	'rmethod','longest');

% Dimension
fem.dim={{'u'}};

% Boundary conditions
fem.border=0;

% Usage
fem.usage=1;

% Problem form
fem.form='coefficient';

% Differentiation
fem.diff='off';

% Differentiation simplification
fem.simplify='on';

% Differentiation rules
fem.rules={};

% Define variables
fem.variables={};

% Define application mode variables
fem.var={};

% Boundary conditions
clear bnd
bnd.q={{{'0'}}};
bnd.g={{{'0'}}};
bnd.h={{{'1'}}};
bnd.r={{{'((x-0.5).^2+(y-0.5).^2).^0.5-0.24'}}};
bnd.var={};
bnd.ind=ones(1,6);
fem.bnd=bnd;

% PDE coefficients
clear equ
equ.da={{{'1./(ux.^2+uy.^2+uz.^2).^0.5'}}};
equ.c={{{'1./(ux.^2+uy.^2+uz.^2).^0.5'}}};
equ.al={{{'0';'0';'0'}}};
equ.ga={{{'0';'0';'0'}}};
equ.be={{{'0';'0';'0'}}};
equ.a={{{'0'}}};
equ.f={{{'0'}}};
equ.var={'absux','sqrt(ux.^2+uy.^2+uz.^2)','abscu1x', ...
'sqrt(cu1x.^2+cu1y.^2+cu1z.^2)'};
equ.ind=1;
fem.equ=equ;

% Evaluate initial condition
fem.init=asseminit(fem,...
	'context','local',...
	'init',   struct('sd',{{{{'((x-0.5).^2+(y-0.5).^2).^0.5-0.24'}}}},'ind',{1}));

% Solve dynamic problem
fem.sol=femtime(fem,...
	'tlist',  0:0.1:1,...
	'atol',   0.001,...
	'rtol',   0.01,...
	'jacobian','equ',...
	'mass',   'full',...
	'ode',    'ode15s',...
	'odeopt', struct('InitialStep',{[]},'MaxOrder',{5},'MaxStep',{[]}),...
	'out',    'sol',...
	'stop',   'on',...
	'report', 'on',...
	'context','local',...
	'sd',     'off',...
	'Epoint', 'gauss2',...
	'Nullfun','flnullorth',...
	'Tpoint', 'gauss2',...
	'Solcomp',1);


% Plot solution
postplot(fem,...
	'context','local',...
	'isocolordata',{'uz','cont','on'},...
	'isodata',{'u','cont','on'},...
	'isofacestyle','flat',...
	'isoedgestyle','none',...
	'isolevels',0,...
	'isostyle','color',...
	'isomaxmin','off',...
	'isobar', 'on',...
	'isomap', 'jet',...
	'geom',   'on',...
	'geomcol','bginv',...
	'lightmodel','phong',...
	'lightreflection','default',...
	'view',   [-10.06577050680178 14.025351442300167],...
	'title',  'Minimal Surface ',...
	'solnum', 11,...
	'renderer','zbuffer',...
	'scenelight','off',...
	'camlight','off',...
	'campos', [-0.7109587412545384 -6.3218868168311815 1.7807337130027319],...
	'camprojection','orthographic',...
	'camtarget',[0.5 0.5 0.05000000000000001],...
	'camup',  [0.037095521954494774 0.23955146431573984 0.97017473590849013],...
	'camva',  7.8620708375802346,...
	'axispos',[-0.18735511064278185 -0.19085173501577285 1.3114857744994732 1.33596214511041],...
	'axisequal','on',...
	'axis',   [-0.050000000000000003 1 -5.5511150999999999e-17 1 -0.050000000000000003 0.20000000000000001],...
	'axisvisible','on',...
	'grid',   'on',...
	'scenelightpos',[-0.7109587412545384 -6.3218868168311815 1.7807337130027319])