rbmatlab  1.16.09
 All Classes Namespaces Files Functions Variables Modules Pages
ThermalBlock_3D_5x5x5_LinLagr_Einfluss_unten_125Params.m
1 function out = model
2 %
3 % ThermalBlock_3D_5x5x5_LinLagr_Einfluss_unten_125Params.m
4 %
5 
6 
7 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
8 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
9 % Set Properties of the model here!!!
10 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
11 %mu_init = [1,1,1]; % [1,3,8,2] --> sets parameters: mu1=1, mu2=3, mu3=8, mu4=2
12 mu_init = ones(125,1);
13 
14 nr_blocks = 125;
15 
16 %zufallszahlen fuer das modell generieren! verwende seed, damit immer die
17 %gleichen zahlen beim model generieren verwendet werden
18 %seed = 1;
19 %rng(seed);
20 %block_mu_index = randi([1 length(mu_init)],nr_blocks,1); %random
21 block_mu_index = 1:length(mu_init);
22 
23 
24 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
25 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
26 % Generation of the Comsol Model
27 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
28 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
29 
30 import com.comsol.model.*
31 import com.comsol.model.util.*
32 
33 model = ModelUtil.create('Model');
34 
35 model.modelPath('/Users/oliverzeeb/Documents/MATLAB/rbmatlab_COMSOL/COMSOL_rbmatlab/models/ThermalBlock/3D_stat');
36 
37 model.name('ThermalBlock_3D_5x5x5_LinLagr_Einfluss_unten_125Params.mph');
38 
39 model.baseSystem('none');
40 
41 for k=1:length(mu_init)
42  model.param.set(['mu',num2str(k)], num2str(mu_init(k)));
43 end
44 
45 model.modelNode.create('mod1');
46 
47 model.geom.create('geom1', 3);
48 %Bloecke generieren
49 for k=1:nr_blocks
50  eval_string = ['model.geom(''geom1'').feature.create(''blk' num2str(k), ''', ''Block'');'];
51  eval(eval_string);
52 end
53 % model.geom('geom1').feature.create('blk1', 'Block');
54 % model.geom('geom1').feature.create('blk2', 'Block');
55 % ...
56 
57 %groesse der bloecke setzen
58 counter = zeros(3,1);
59 for k=1:nr_blocks
60  %position festlegen
61  pos_string = ['model.geom(''geom1'').feature(''blk', num2str(k) , ''').set(''pos'', {'''...
62  num2str(counter(1)) '/5'' ''' num2str(counter(2)), '/5'' ''', num2str(counter(3)) '/5''});'];
63  eval(pos_string);
64 
65  size_string = ['model.geom(''geom1'').feature(''blk', num2str(k), ''').set(''size'', {''1/5'' ''1/5'' ''1/5''});'];
66  eval(size_string);
67 
68  if k<nr_blocks
69  counter = increase_counter(5,1,counter);
70  end
71 end
72 % model.geom('geom1').feature('blk1').set('pos', {'0' '0' '0'});
73 % model.geom('geom1').feature('blk1').set('size', {'1/4' '1/4' '1/4'});
74 % ...
75 
76 model.geom('geom1').run;
77 
78 model.physics.create('c', 'CoefficientFormPDE', 'geom1');
79 for k=2:nr_blocks
80  create_string = ['model.physics(''c'').feature.create(''cfeq' num2str(k) ''', ''CoefficientFormPDE'', 3);'];
81  eval(create_string);
82  selection_string = ['model.physics(''c'').feature(''cfeq' num2str(k) ''').selection.set([' num2str(k) ']);'];
83  eval(selection_string);
84 end
85 % model.physics('c').feature.create('cfeq2', 'CoefficientFormPDE', 3);
86 % model.physics('c').feature('cfeq2').selection.set([2]);
87 % ...
88 
89 model.physics('c').feature.create('dir1', 'DirichletBoundary', 2);
90 model.physics('c').feature('dir1').selection.set([16, 32, 48, 64, 80, 101, 117, 133, 149, 165,...
91  186, 202, 218, 234, 250, 271, 287, 303, 319, 335, 356, 372, 388, 404, 420]);
92 model.physics('c').feature.create('flux1', 'FluxBoundary', 2);
93 model.physics('c').feature('flux1').selection.set([3, 19, 35, 51, 67, 88, 104, 120, 136, 152,...
94  173, 189, 205, 221, 237, 258, 274, 290, 306, 322, 343, 359, 375, 391, 407]);
95 
96 model.mesh.create('mesh1', 'geom1');
97 model.mesh('mesh1').feature.create('ftet1', 'FreeTet');
98 
99 model.view('view1').set('renderwireframe', true);
100 
101 model.physics('c').prop('ShapeProperty').set('order', '1');
102 
103 % koeffizienten der PDE setzen
104 for k=1:nr_blocks
105  set_f_string = ['model.physics(''c'').feature(''cfeq' num2str(k) ''').set(''f'', ''0'');'];
106  eval(set_f_string)
107  set_da_string = ['model.physics(''c'').feature(''cfeq' num2str(k) ''').set(''da'', ''0'');'];
108  eval(set_da_string)
109 end
110 model.physics('c').feature('flux1').set('g', '50');
111 for k=1:length(block_mu_index)
112  mu_string = ['mu',num2str(block_mu_index(k))];
113  model.physics('c').feature(['cfeq',num2str(k)]).set('c', {mu_string '0' '0' '0' mu_string '0' '0' '0' mu_string});
114 end
115 % model.physics('c').feature('cfeq1').set('c', {'mu1' '0' '0' '0' 'mu1' '0' '0' '0' 'mu1'});
116 % model.physics('c').feature('cfeq1').set('f', '0');
117 % model.physics('c').feature('cfeq1').set('da', '0');
118 % ...
119 
120 %model.mesh('mesh1').run;
121 
122 model.frame('material1').sorder(1);
123 
124 model.study.create('std1');
125 model.study('std1').feature.create('stat', 'Stationary');
126 
127 model.sol.create('sol1');
128 model.sol('sol1').study('std1');
129 model.sol('sol1').attach('std1');
130 model.sol('sol1').feature.create('st1', 'StudyStep');
131 model.sol('sol1').feature.create('v1', 'Variables');
132 model.sol('sol1').feature.create('s1', 'Stationary');
133 model.sol('sol1').feature('s1').feature.create('fc1', 'FullyCoupled');
134 model.sol('sol1').feature('s1').feature.remove('fcDef');
135 
136 % Ausgabefuntional
137 model.result.numerical.create('int1', 'IntSurface');
138 model.result.numerical('int1').selection.set([3, 19, 35, 51, 67, 88, 104, 120, 136, 152,...
139  173, 189, 205, 221, 237, 258, 274, 290, 306, 322, 343, 359, 375, 391, 407]);
140 model.result.numerical('int1').set('probetag', 'none');
141 
142 model.result.create('pg1', 'PlotGroup3D');
143 model.result('pg1').feature.create('slc1', 'Slice');
144 
145 model.sol('sol1').feature('st1').name('Compile Equations: Stationary {stat}');
146 model.sol('sol1').feature('st1').set('studystep', 'stat');
147 model.sol('sol1').feature('v1').set('control', 'stat');
148 model.sol('sol1').feature('s1').set('control', 'stat');
149 %model.sol('sol1').runAll;
150 
151 model.result('pg1').feature('slc1').set('unit', '1');
152 
153 out = model;
154 
155 end
156 
157 
158 
159 
160 
161 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
162 % Auxiliary functions
163 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
164 function counter = increase_counter(basis, pos,counter)
165  if mod(counter(pos)+1,basis) == 0 %aktuelle pos auf null setzen, nächste pos um 1 erhöhen!
166  counter(pos) = 0;
167  counter = increase_counter(basis, pos+1,counter);
168  else
169  counter(pos) = counter(pos)+1;
170  end
171 end