function out = model
%
% ModBenchmark.m
%
% Model exported on Jan 28 2014, 15:06 by COMSOL 4.3.2.189.

import com.comsol.model.*
import com.comsol.model.util.*

model = ModelUtil.create('Model');

%model.modelPath('/home/s1061559/MEngProject');

% mod1 = model.modelNode.create('mod1');

geom1 = model.geom.create('geom1', 2);

mesh1 = model.mesh.create('mesh1', 'geom1');

slpf = model.physics.create('slpf', 'SlipFlow', 'geom1');

ModelUtil.showProgress(true)

% model.study.create('std1');
% model.study('std1').feature.create('stat', 'Stationary');
% model.study('std1').feature('stat').activate('slpf', true);
std1 = model.study.create('std1');
stat = std1.feature.create('stat', 'Stationary');
stat.activate('slpf', true);

% ***** Create Global Parameters *****
% Geometry Parameters
model.param.set('b', '10[um]', 'Length of the reservoir');
model.param.set('d', '6.25[um]', 'Depth of the reservoir');
model.param.set('D', '1.5[um]', 'Channel width');
model.param.set('Lc', '15[um]', 'Channel length');
model.param.set('Lp', 'Lc+b', 'Pump length');
model.param.set('Rh', '2*d+D', 'Reservoir height');
model.param.set('tw', '1[um]', 'Wall thickness');

% Property Parameters
model.param.set('Tc', '300[K]', 'Temperature cold');
model.param.set('Th', '800[K]', 'Temperature hot');

% ***** Create The Geometry *****
% Create Geometry Node

% Create the Shell
r1 = geom1.feature.create('r1', 'Rectangle');
r1.set('size', {'(Lp-Lc)/2+tw' 'Rh+2*tw'});
r2 = geom1.feature.create('r2', 'Rectangle');
r2.set('size', {'Lc' 'D+2*tw'});
r2.set('pos', {'b/2' 'd'});
r3 = geom1.feature.create('r3', 'Rectangle');
r3.set('size', {'(Lp-Lc)/2+tw' 'Rh+2*tw'});
r3.set('pos', {'b/2+Lc-tw' '0'});
% Create geomety unions
uni1 = geom1.feature.create('uni1', 'Union');
uni1.selection('input').set({'r1' 'r2' 'r3'});
uni1.set('intbnd', 'off');

% Create the Interior
r4 = geom1.feature.create('r4', 'Rectangle');
r4.set('size', {'b/2' 'Rh'});
r4.set('pos', {'0' 'tw'});
r5 = geom1.feature.create('r5', 'Rectangle');
r5.set('size', {'Lc' 'D'});
r5.set('pos', {'b/2' 'd+tw'});
r6 = geom1.feature.create('r6', 'Rectangle');
r6.set('size', {'b/2' 'Rh'});
r6.set('pos', {'b/2+Lc' 'tw'});
% Unite the interior geometry
uni2 = geom1.feature.create('uni2', 'Union');
uni2.selection('input').set({'r4' 'r5' 'r6'});
uni2.set('intbnd', 'off');
mat = model.material.create('mat_air');

% Create pressure point constraint
p1 = geom1.feature.create('p1', 'Point');
p1.set('p', {'Lp/2' '(Rh+2*tw)/2'});
geom1.run

% Create a pair
pr1 = model.pair.create('pr1', 'Identity', 'geom1');
pr1.source.set([3 25]);
pr1.destination.set([3 25]);

geom1.run
mphgeom(model, 'geom1', 'Facelabels', 'on', 'Edgelabels', 'on')

%% Old Geometry
% model.geom('geom1').lengthUnit('km');
% model.geom('geom1').lengthUnit([native2unicode(hex2dec('00b5'), 'Cp1252') 'm']);
% model.geom('geom1').feature.create('r1', 'Rectangle');
% model.geom('geom1').feature('r1').setIndex('size', '15', 0);
% model.geom('geom1').feature('r1').setIndex('size', '1.5', 1);
% model.geom('geom1').feature('r1').setIndex('pos', '-.75', 1);
% model.geom('geom1').run('r1');
% model.geom('geom1').run('r1');
% model.geom('geom1').feature.create('r2', 'Rectangle');
% model.geom('geom1').feature('r2').setIndex('size', '7.5', 0);
% model.geom('geom1').feature('r2').setIndex('size', '15', 1);
% model.geom('geom1').feature('r2').setIndex('pos', '-7.5', 0);
% model.geom('geom1').feature('r2').setIndex('pos', '-7.5', 1);
% model.geom('geom1').run('r2');
% model.geom('geom1').run('r2');
% model.geom('geom1').feature.create('r3', 'Rectangle');
% model.geom('geom1').feature('r3').setIndex('size', '7.5', 0);
% model.geom('geom1').feature('r3').setIndex('size', '15', 1);
% model.geom('geom1').feature('r3').setIndex('pos', '15', 0);
% model.geom('geom1').feature('r3').setIndex('pos', '-7.5', 1);
% model.geom('geom1').run('r3');
% model.geom('geom1').run('r3');
% model.geom('geom1').feature.create('uni1', 'Union');
% model.geom('geom1').feature('uni1').selection('input').set({'r1' 'r2' 'r3'});
% model.geom('geom1').runAll;
% model.geom('geom1').feature('uni1').set('intbnd', 'off');
% model.geom('geom1').runAll;
% model.geom('geom1').run('uni1');
% model.geom('geom1').feature.create('r4', 'Rectangle');
% model.geom('geom1').feature('r4').setIndex('size', '8.6', 0);
% model.geom('geom1').feature('r4').setIndex('size', '17', 1);
% model.geom('geom1').feature('r4').setIndex('pos', '-7.5', 0);
% model.geom('geom1').feature('r4').setIndex('pos', '-8.5', 1);
% model.geom('geom1').run('r4');
% model.geom('geom1').run('r4');
% model.geom('geom1').feature.create('r5', 'Rectangle');
% model.geom('geom1').feature('r5').setIndex('size', '15', 0);
% model.geom('geom1').feature('r5').setIndex('size', '3.5', 1);
% model.geom('geom1').feature('r5').setIndex('pos', '-1.75', 1);
% model.geom('geom1').run('r5');
% model.geom('geom1').run('r5');
% model.geom('geom1').feature.create('r6', 'Rectangle');
% model.geom('geom1').feature('r6').setIndex('size', '8.5', 0);
% model.geom('geom1').feature('r6').setIndex('size', '17', 1);
% model.geom('geom1').feature('r6').setIndex('pos', '14', 0);
% model.geom('geom1').feature('r6').setIndex('pos', '-8.5', 1);
% model.geom('geom1').run('r6');
% model.geom('geom1').run('r6');
% model.geom('geom1').feature.create('uni2', 'Union');
% model.geom('geom1').feature('uni2').selection('input').set({'r4' 'r5' 'r6'});
% model.geom('geom1').run('uni2');
% model.geom('geom1').feature('uni2').set('intbnd', 'off');
% model.geom('geom1').run('uni2');
% model.geom('geom1').run('uni2');
% model.geom('geom1').feature.create('pt1', 'Point');
% model.geom('geom1').feature('pt1').setIndex('p', '7.5', 0);
% model.geom('geom1').run('pt1');
% model.geom('geom1').run;

%% Material Details
% Material 1 - Air
mat1 = model.material.create('mat1');
mat1.set('family', 'air');
mat1.materialModel('def').set('relpermeability', '1');
mat1.materialModel('def').set('relpermittivity', '1');
mat1.materialModel('def').set('dynamicviscosity', 'eta(T[1/K])[Pa*s]');
mat1.materialModel('def').set('ratioofspecificheat', '1.4');
mat1.materialModel('def').set('electricconductivity', '0[S/m]');
mat1.materialModel('def').set('heatcapacity', 'Cp(T[1/K])[J/(kg*K)]');
mat1.materialModel('def').set('density', 'rho(pA[1/Pa],T[1/K])[kg/m^3]');
mat1.materialModel('def').set('thermalconductivity', 'k(T[1/K])[W/(m*K)]');
mat1.materialModel('def').set('soundspeed', 'cs(T[1/K])[m/s]');
mat1.materialModel('def').func.create('eta', 'Piecewise');
mat1.materialModel('def').func('eta').set('funcname', 'eta');
mat1.materialModel('def').func('eta').set('arg', 'T');
mat1.materialModel('def').func('eta').set('extrap', 'constant');
mat1.materialModel('def').func('eta').set('pieces', {'200.0' '1600.0' '-8.38278E-7+8.35717342E-8*T^1-7.69429583E-11*T^2+4.6437266E-14*T^3-1.06585607E-17*T^4'});
mat1.materialModel('def').func.create('Cp', 'Piecewise');
mat1.materialModel('def').func('Cp').set('funcname', 'Cp');
mat1.materialModel('def').func('Cp').set('arg', 'T');
mat1.materialModel('def').func('Cp').set('extrap', 'constant');
mat1.materialModel('def').func('Cp').set('pieces', {'200.0' '1600.0' '1047.63657-0.372589265*T^1+9.45304214E-4*T^2-6.02409443E-7*T^3+1.2858961E-10*T^4'});
mat1.materialModel('def').func.create('rho', 'Analytic');
mat1.materialModel('def').func('rho').set('funcname', 'rho');
mat1.materialModel('def').func('rho').set('args', {'pA' 'T'});
mat1.materialModel('def').func('rho').set('expr', 'pA*0.02897/8.314/T');
mat1.materialModel('def').func('rho').set('dermethod', 'manual');
mat1.materialModel('def').func('rho').set('argders', {'pA' 'd(pA*0.02897/8.314/T,pA)'; 'T' 'd(pA*0.02897/8.314/T,T)'});
mat1.materialModel('def').func.create('k', 'Piecewise');
mat1.materialModel('def').func('k').set('funcname', 'k');
mat1.materialModel('def').func('k').set('arg', 'T');
mat1.materialModel('def').func('k').set('extrap', 'constant');
mat1.materialModel('def').func('k').set('pieces', {'200.0' '1600.0' '-0.00227583562+1.15480022E-4*T^1-7.90252856E-8*T^2+4.11702505E-11*T^3-7.43864331E-15*T^4'});
mat1.materialModel('def').func.create('cs', 'Analytic');
mat1.materialModel('def').func('cs').set('funcname', 'cs');
mat1.materialModel('def').func('cs').set('args', {'T'});
mat1.materialModel('def').func('cs').set('expr', 'sqrt(1.4*287*T)');
mat1.materialModel('def').func('cs').set('dermethod', 'manual');
mat1.materialModel('def').func('cs').set('argders', {'T' 'd(sqrt(1.4*287*T),T)'});
mat1.materialModel('def').addInput('temperature');
mat1.materialModel('def').addInput('pressure');
mat1.materialModel.create('RefractiveIndex', 'Refractive index');
mat1.materialModel('RefractiveIndex').set('n', '1');
mat1.materialModel('def').set('family', 'air');
mat1.name('Air');

%% Material 2 - Silicon
mat2 = model.material.create('mat2');
mat2.name('Silicon');
mat2.set('family', 'custom');

mat2.set('lighting', 'cooktorrance');
mat2.set('specular', 'custom');
mat2.set('customspecular', [0.7843137254901961 1 1]);

mat2.set('fresnel', 0.9);
mat2.set('roughness', 0.1);
mat2.set('shininess', 200);
mat2.set('diffuse', 'custom');
mat2.set('customdiffuse', [0.6666666666666666 0.6666666666666666 0.7058823529411765]);
mat2.set('ambient', 'custom');
mat2.set('customambient', [0.6666666666666666 0.6666666666666666 0.7058823529411765]);
mat2.set('fresnel', 0.7);
mat2.set('roughness', 0.5);
mat2.materialModel('def').set('relpermeability', '1');
mat2.materialModel('def').set('electricconductivity', '1e-12[S/m]');
mat2.materialModel('def').set('thermalexpansioncoefficient', '2.6e-6[1/K]');
mat2.materialModel('def').set('heatcapacity', '700[J/(kg*K)]');
mat2.materialModel('def').set('relpermittivity', '11.7');
mat2.materialModel('def').set('density', '2329[kg/m^3]');
mat2.materialModel('def').set('thermalconductivity', '130[W/(m*K)]');
mat2.materialModel.create('Enu', 'Young''s modulus and Poisson''s ratio');
mat2.materialModel('Enu').set('poissonsratio', '0.28');
mat2.materialModel('Enu').set('youngsmodulus', '170e9[Pa]');
mat2.materialModel.create('RefractiveIndex', 'Refractive index');
mat2.materialModel('RefractiveIndex').set('n', '3.48');
mat2.selection.set([1 3]);

%% Old Mat1
% model.material.create('mat1');
% model.material('mat1').name('Air');
% model.material('mat1').set('family', 'air');
% model.material('mat1').propertyGroup('def').set('relpermeability', '1');
% model.material('mat1').propertyGroup('def').set('relpermittivity', '1');
% model.material('mat1').propertyGroup('def').set('dynamicviscosity', 'eta(T[1/K])[Pa*s]');
% model.material('mat1').propertyGroup('def').set('ratioofspecificheat', '1.4');
% model.material('mat1').propertyGroup('def').set('electricconductivity', '0[S/m]');
% model.material('mat1').propertyGroup('def').set('heatcapacity', 'Cp(T[1/K])[J/(kg*K)]');
% model.material('mat1').propertyGroup('def').set('density', 'rho(pA[1/Pa],T[1/K])[kg/m^3]');
% model.material('mat1').propertyGroup('def').set('thermalconductivity', 'k(T[1/K])[W/(m*K)]');
% model.material('mat1').propertyGroup('def').set('soundspeed', 'cs(T[1/K])[m/s]');
% model.material('mat1').propertyGroup('def').func.create('eta', 'Piecewise');
% model.material('mat1').propertyGroup('def').func('eta').set('funcname', 'eta');
% model.material('mat1').propertyGroup('def').func('eta').set('arg', 'T');
% model.material('mat1').propertyGroup('def').func('eta').set('extrap', 'constant');
% model.material('mat1').propertyGroup('def').func('eta').set('pieces', {'200.0' '1600.0' '-8.38278E-7+8.35717342E-8*T^1-7.69429583E-11*T^2+4.6437266E-14*T^3-1.06585607E-17*T^4'});
% model.material('mat1').propertyGroup('def').func.create('Cp', 'Piecewise');
% model.material('mat1').propertyGroup('def').func('Cp').set('funcname', 'Cp');
% model.material('mat1').propertyGroup('def').func('Cp').set('arg', 'T');
% model.material('mat1').propertyGroup('def').func('Cp').set('extrap', 'constant');
% model.material('mat1').propertyGroup('def').func('Cp').set('pieces', {'200.0' '1600.0' '1047.63657-0.372589265*T^1+9.45304214E-4*T^2-6.02409443E-7*T^3+1.2858961E-10*T^4'});
% model.material('mat1').propertyGroup('def').func.create('rho', 'Analytic');
% model.material('mat1').propertyGroup('def').func('rho').set('funcname', 'rho');
% model.material('mat1').propertyGroup('def').func('rho').set('args', {'pA' 'T'});
% model.material('mat1').propertyGroup('def').func('rho').set('expr', 'pA*0.02897/8.314/T');
% model.material('mat1').propertyGroup('def').func('rho').set('dermethod', 'manual');
% model.material('mat1').propertyGroup('def').func('rho').set('argders', {'pA' 'd(pA*0.02897/8.314/T,pA)'; 'T' 'd(pA*0.02897/8.314/T,T)'});
% model.material('mat1').propertyGroup('def').func.create('k', 'Piecewise');
% model.material('mat1').propertyGroup('def').func('k').set('funcname', 'k');
% model.material('mat1').propertyGroup('def').func('k').set('arg', 'T');
% model.material('mat1').propertyGroup('def').func('k').set('extrap', 'constant');
% model.material('mat1').propertyGroup('def').func('k').set('pieces', {'200.0' '1600.0' '-0.00227583562+1.15480022E-4*T^1-7.90252856E-8*T^2+4.11702505E-11*T^3-7.43864331E-15*T^4'});
% model.material('mat1').propertyGroup('def').func.create('cs', 'Analytic');
% model.material('mat1').propertyGroup('def').func('cs').set('funcname', 'cs');
% model.material('mat1').propertyGroup('def').func('cs').set('args', {'T'});
% model.material('mat1').propertyGroup('def').func('cs').set('expr', 'sqrt(1.4*287*T)');
% model.material('mat1').propertyGroup('def').func('cs').set('dermethod', 'manual');
% model.material('mat1').propertyGroup('def').func('cs').set('argders', {'T' 'd(sqrt(1.4*287*T),T)'});
% model.material('mat1').propertyGroup('def').addInput('temperature');
% model.material('mat1').propertyGroup('def').addInput('pressure');
% model.material('mat1').propertyGroup.create('RefractiveIndex', 'Refractive index');
% model.material('mat1').propertyGroup('RefractiveIndex').set('n', '1');
% model.material('mat1').set('family', 'air');

%% Old Mat2
% model.material.create('mat2');
% model.material('mat2').name('Silicon');
% model.material('mat2').set('family', 'custom');
% model.material('mat2').set('lighting', 'cooktorrance');
% model.material('mat2').set('specular', 'custom');
% model.material('mat2').set('customspecular', [0.7843137254901961 1 1]);
% model.material('mat2').set('fresnel', 0.9);
% model.material('mat2').set('roughness', 0.1);
% model.material('mat2').set('shininess', 200);
% model.material('mat2').set('diffuse', 'custom');
% model.material('mat2').set('customdiffuse', [0.6666666666666666 0.6666666666666666 0.7058823529411765]);
% model.material('mat2').set('ambient', 'custom');
% model.material('mat2').set('customambient', [0.6666666666666666 0.6666666666666666 0.7058823529411765]);
% model.material('mat2').set('fresnel', 0.7);
% model.material('mat2').set('roughness', 0.5);
% model.material('mat2').propertyGroup('def').set('relpermeability', '1');
% model.material('mat2').propertyGroup('def').set('electricconductivity', '1e-12[S/m]');
% model.material('mat2').propertyGroup('def').set('thermalexpansioncoefficient', '2.6e-6[1/K]');
% model.material('mat2').propertyGroup('def').set('heatcapacity', '700[J/(kg*K)]');
% model.material('mat2').propertyGroup('def').set('relpermittivity', '11.7');
% model.material('mat2').propertyGroup('def').set('density', '2329[kg/m^3]');
% model.material('mat2').propertyGroup('def').set('thermalconductivity', '130[W/(m*K)]');
% model.material('mat2').propertyGroup.create('Enu', 'Young''s modulus and Poisson''s ratio');
% model.material('mat2').propertyGroup('Enu').set('poissonsratio', '0.28');
% model.material('mat2').propertyGroup('Enu').set('youngsmodulus', '170e9[Pa]');
% model.material('mat2').propertyGroup.create('RefractiveIndex', 'Refractive index');
% model.material('mat2').propertyGroup('RefractiveIndex').set('n', '3.48');
% model.material('mat2').set('family', 'custom');
% model.material('mat2').set('lighting', 'cooktorrance');
% model.material('mat2').set('specular', 'custom');
% model.material('mat2').set('customspecular', [0.7843137254901961 1 1]);
% model.material('mat2').set('fresnel', 0.9);
% model.material('mat2').set('roughness', 0.1);
% model.material('mat2').set('shininess', 200);
% model.material('mat2').set('diffuse', 'custom');
% model.material('mat2').set('customdiffuse', [0.6666666666666666 0.6666666666666666 0.7058823529411765]);
% model.material('mat2').set('ambient', 'custom');
% model.material('mat2').set('customambient', [0.6666666666666666 0.6666666666666666 0.7058823529411765]);
% model.material('mat2').set('fresnel', 0.7);
% model.material('mat2').set('roughness', 0.5);
% model.material('mat2').selection.set([1 3]);

%% Physics
solid1 = slpf.feature.create('solid1', 'SolidHeatTransferModel', 2);
solid1.selection.set([1 3])
ins2 = slpf.feature.create('ins2', 'ThermalInsulation', 1);
ins2.selection.set([13 15]);
temp1 = slpf.feature.create('temp1', 'TemperatureBoundary', 1);
temp1.selection.set([1 2 5 7 12 14])
temp1.set('T0', 1, 'Tc')
temp2 = slpf.feature.create('temp2', 'TemperatureBoundary', 1);
temp2.selection.set([16:19 24 26])
temp2.set('T0', 1, 'Th')
slw1 = slpf.feature.create('slw1', 'SlipWall', 1);
slw1.selection.set([8:11 20 22])
prpc1 = slpf.feature.create('prpc1', 'PressurePointConstraint', 0);
prpc1.selection.set([13]);
pfc1 = slpf.feature.create('pfc1', 'PeriodicFlowCondition', 1);
pfc1.selection.set([3 25]);
tempb1 = slpf.feature.create('tempb1', 'PairTemperatureBoundary', 1);
tempb1.set('T0', 1, '550[K]');
tempb1.set('pairs', 1, 'pr1');
%% Old Physics
% model.physics('slpf').feature.create('solid1', 'SolidHeatTransferModel', 2);
% model.physics('slpf').feature('solid1').selection.set([1 3]);
% model.physics('slpf').feature.create('ins2', 'ThermalInsulation', 1);
% model.physics('slpf').feature('ins2').selection.set([13 15]);
% model.physics('slpf').feature.create('temp1', 'TemperatureBoundary', 1);
% model.physics('slpf').feature('temp1').selection.set([1 2 5 7 12 14]);
% model.physics('slpf').feature('temp1').set('T0', 1, '300[K]');
% model.physics('slpf').feature.create('temp2', 'TemperatureBoundary', 1);
% model.physics('slpf').feature('temp2').set('T0', 1, '800[K]');
% model.physics('slpf').feature('temp2').selection.set([16 17 18 19 24 26]);
% model.physics('slpf').feature.create('slw1', 'SlipWall', 1);
% model.physics('slpf').feature('slw1').selection.set([8 9 10 11 20 22]);
% model.physics('slpf').feature.create('prpc1', 'PressurePointConstraint', 0);
% model.physics('slpf').feature('prpc1').selection.set([13]);
% 
% model.pair.create('p1', 'Identity', 'geom1');
% model.pair('p1').source.set([3 25]);
% model.pair('p1').destination.set([3 25]);
% 
% model.physics('slpf').feature.create('temp3', 'PairTemperatureBoundary', 1);
% model.physics('slpf').feature('temp3').set('pairs', 1, 'p1');
% model.physics('slpf').feature('temp3').set('T0', 1, '550[K]');
% model.physics('slpf').feature.create('pfc1', 'PeriodicFlowCondition', 1);
% model.physics('slpf').feature('pfc1').selection.set([3 25]);

%% Meshing

size = mesh1.feature('size')
size.set('table', 'cfd');
size.set('hauto', 1);
ftri1 = mesh1.feature.create('ftri1', 'FreeTri');
mesh1.run
% model.mesh('mesh1').feature.create('size1', 'Size');
% model.mesh('mesh1').feature('size1').set('table', 'cfd');
% model.mesh('mesh1').feature('size1').set('hauto', '2');
% model.mesh('mesh1').feature.create('ftri1', 'FreeTri');
% model.mesh('mesh1').run;
mphgeom(model)

std1.run

%% Results
% model.result.create('pg1', 2);
% model.result('pg1').set('data', 'dset1');
% model.result('pg1').feature.create('surf1', 'Surface');
% model.result('pg1').feature('surf1').set('expr', {'slpf.U'});
% model.result('pg1').set('frametype', 'spatial');
% model.result('pg1').name('Velocity (slpf)');
% model.result.create('pg2', 2);
% model.result('pg2').set('data', 'dset1');
% model.result('pg2').feature.create('con1', 'Contour');
% model.result('pg2').feature('con1').set('expr', {'p'});
% model.result('pg2').set('frametype', 'spatial');
% model.result('pg2').name('Pressure (slpf)');
% model.result('pg2').feature('con1').set('number', 40);
% 
% model.sol('sol1').runAll;
% 
% model.result('pg1').run;
% model.result.create('pg3', 'PlotGroup2D');
% model.result('pg3').run;
% model.result('pg3').name('Mean Free Path');
% model.result('pg3').feature.create('surf1', 'Surface');
% model.result('pg3').feature('surf1').set('expr', 'slpf.lambda');
% model.result('pg3').run;
% 
% 
% 
% 
% model.result('pg1').run;
% mphplot(model, 'pg1', 'rangenum', 1)

out = model;