% COMSOL Multiphysics Model M-file
% Generated by COMSOL 3.2b (COMSOL 3.2.0.304, $Date: 2006/04/04 14:56:13 $)
function [arrayout] = FieldFile (varargin)
flclear fem

% COMSOL version
clear vrsn
vrsn.name = 'COMSOL 3.2';
vrsn.ext = 'b';
vrsn.major = 0;
vrsn.build = 304;
vrsn.rcs = '$Name:  $';
vrsn.date = '$Date: 2006/04/04 14:56:13 $';
fem.version = vrsn;

% Geometry
g1=rect2(3,750,'base','corner','pos',[-1,0]);
a = length(varargin);
a;
x = varargin{1}; % Cell array indexing
y = varargin{2};
V = varargin(3);
VV = V{1}
x';
y';
% x=[0 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20];
% y=sin(x/20*pi)/100;

carr={curve2([x(1),x(2)],[y(1),y(2)],[1,1]), ...
  curve2([x(2),x(3)],[y(2),y(3)],[1,1]), ...
  curve2([x(3),x(4)],[y(3),y(4)],[1,1]), ...
  curve2([x(4),x(5)],[y(4),y(5)],[1,1]), ...
  curve2([x(5),x(6)],[y(5),y(6)],[1,1]), ...
  curve2([x(6),x(7)],[y(6),y(7)],[1,1]), ...
  curve2([x(7),x(8)],[y(7),y(8)],[1,1]), ...
  curve2([x(8),x(9)],[y(8),y(9)],[1,1]), ...
  curve2([x(9),x(10)],[y(9),y(10)],[1,1]), ...
  curve2([x(10),x(11)],[y(10),y(11)],[1,1]), ...
  curve2([x(11),x(12)],[y(11),y(12)],[1,1]), ...
  curve2([x(12),x(13)],[y(12),y(13)],[1,1]), ...
  curve2([x(13),x(14)],[y(13),y(14)],[1,1]), ...
  curve2([x(14),x(15)],[y(14),y(15)],[1,1]), ...
  curve2([x(15),x(16)],[y(15),y(16)],[1,1]), ...
  curve2([x(16),x(17)],[y(16),y(17)],[1,1]), ...
  curve2([x(17),x(18)],[y(17),y(18)],[1,1]), ...
  curve2([x(18),x(19)],[y(18),y(19)],[1,1]), ...
  curve2([x(19),x(20)],[y(19),y(20)],[1,1]), ...
  curve2([x(20),x(21)],[y(20),y(21)],[1,1]), ...
  curve2([x(21),x(1)],[y(21),y(1)],[1,1])};
g2=geomcoerce('solid',carr);
g3=geomcomp({g1,g2},'ns',{'g1','g2'},'sf','g1+g2','edge','none');
garr=split(g3);
[g4,g5]=deal(garr{:});
clear s
s.objs={g4,g5};
s.name={'CO1','CO3'};
s.tags={'g4','g5'};

fem.draw=struct('s',s);
fem.geom=geomcsg(fem);

% Initialize mesh
fem.mesh=meshinit(fem);

% (Default values are not included)

% Application mode 1
clear appl
appl.mode.class = 'EmElectrostatics';
appl.module = 'EM';
appl.assignsuffix = '_emes';
clear bnd
bnd.V0 = {0,0,-VV,VV};
bnd.type = {'V0','cont','V','V'};
bnd.ind = [1,3,4,2,3,2,2,2,2,2,2,2,2,2,2,2,2,2,2,2,2,2,2,2,3,1];
appl.bnd = bnd;
clear equ
equ.epsilonr = {1e9,5e9};
equ.ind = [1,2];
appl.equ = equ;
fem.appl{1} = appl;
fem.frame = {'ref'};
fem.border = 1;
fem.units = 'SI';

% Multiphysics
fem=multiphysics(fem);

% Extend mesh
fem.xmesh=meshextend(fem);

% Solve problem
fem.sol=femlin(fem, ...
               'solcomp',{'V'}, ...
               'outcomp',{'V'});

% Save current fem structure for restart purposes
%fem0=fem;

%Plot solution

% postplot(fem, ...
%          'tridata',{'We_emes','cont','internal'}, ...
%          'trimap','jet(1024)', ...
%          'title','Surface: Component [We_emes]', ...
%          'axis',[-2.50480121262003,5.49519921262003,-1.053211,2.297248,-1,1]);

 gridd1 = [x(:),y(:)+0.5*y(:)]';
 
 mat_result1 = postinterp(fem,'We_emes',gridd1);
 w = mat_result1;
 

arrayout = {w};
flclear fem;
end