MoDeNa  1.0
Software framework facilitating sequential multi-scale modelling
MomSources.H
Go to the documentation of this file.
1 
20 if (PBESwitch)
22 {
23  if
24  (
25  PBEMethod == "QMOM"
26  && nNodes == 2
27  )
28  {
29  static int nodes = nNodes;
30  int mOrder[2*nodes];
31  for (int i = 0; i < 2*nNodes; i++){mOrder[i] = i;}
32  double mom[2*nodes], we[nodes], vi[nodes],
33  sgBA[2*nodes], sgCO2[2*nodes], sc[2*nodes];
34  double sources[2*nodes];
35 
36  double radiusGrowthBA[nodes],
37  radiusGrowthCO2[nodes],
38  volumeGrowthBA[nodes],
39  volumeGrowthCO2[nodes],
40  nodeRadii[nodes];
41  int ret_bblgr1[nodes], ret_bblgr2[nodes];
42 
43  double CO2_l_val, L_l_val, tmp_val, alphafoam;
44  double cc1_val = 0.0;
45  double LMAX, xBL_val, rhoPolyS;
46  // Calling the model for density reaction mixture
47  size_t Td_pos =
48  modena_model_inputs_argPos(density_reaction_mixturemodel, "T");
49  size_t XOHd_pos =
50  modena_model_inputs_argPos(density_reaction_mixturemodel, "XOH");
51  modena_model_argPos_check(density_reaction_mixturemodel);
52 
53  // Calling the bubble surrogate models
54  size_t Tbblgr1pos = modena_model_inputs_argPos(bblgr1, "T");
55  size_t Rbblgr1pos = modena_model_inputs_argPos(bblgr1, "R");
56  // size_t KH1bblgr1pos = modena_model_inputs_argPos(bblgr1, "kH");
57  size_t c_1bblgr1pos = modena_model_inputs_argPos(bblgr1, "c");
58  size_t p_1bblgr1pos = modena_model_inputs_argPos(bblgr1, "p");
60 
61  size_t Tbblgr2pos = modena_model_inputs_argPos(bblgr2, "T");
62  size_t Rbblgr2pos = modena_model_inputs_argPos(bblgr2, "R");
63  // size_t KH2bblgr2pos = modena_model_inputs_argPos(bblgr2, "kH");
64  size_t c_2bblgr2pos = modena_model_inputs_argPos(bblgr2, "c");
65  size_t p_2bblgr2pos = modena_model_inputs_argPos(bblgr2, "p");
67 
68  forAll(mesh.C(), celli)
69  {
70  CO2_l_val = wCO2_l[celli];
71  L_l_val = wBA_l[celli];
72  tmp_val = TS[celli];
73  alphafoam = alpha2[celli];
74  double wCO2_g_val, wBA_g_val,prsr;
75  wCO2_g_val= wCO2_g[celli];
76  wBA_g_val = wBA_g[celli];
77  prsr = p[celli];
78 
79  if (liquidMixtureDensitySurrogate)
80  {
81  // set input vector
82  modena_inputs_set(inputs_den, Td_pos, TS[celli]);
83  modena_inputs_set(inputs_den, XOHd_pos, XOH[celli]);
84  int ret_den =
86  (
87  density_reaction_mixturemodel,
88  inputs_den,
89  outputs_den
90  );
91  if (ret_den != 0)
92  {
93  modena_inputs_destroy (inputs_den);
94  modena_outputs_destroy (outputs_den);
95  modena_model_destroy (density_reaction_mixturemodel);
96  exit(ret_den);
97  }
98  rhoPolyS = modena_outputs_get(outputs_den, 0);
99  }
100  else
101  {
102  rhoPolyS = rhoPoly;
103  }
104 
105  // Using old moments to compute the new weights and nodes
106  if (alphafoam > 0.5)
107  {
108  mom[0] = mZero[celli];
109  mom[1] = mOne[celli];
110  mom[2] = mTwo[celli];
111  mom[3] = mThree[celli];
112 
113  // calling the functions, PDA, growthSource, coalescenceSource
114  PDA(we, vi, mom, nodes);
115  weight0[celli] = we[0];
116  weight1[celli] = we[1];
117  node0[celli] = vi[0];
118  node1[celli] = vi[1];
119 
120  // Required variables for wCO2_Max
121  double henry_Coefficient, bubble_radius,
122  partialPressure_CO2, wCO2_Max;
123  henry_Coefficient = henryCoefficient(tmp_val);
124  bubble_radius = bubbleRadius(mom[0], mom[1]);
125  partialPressure_CO2 =
127  (
128  M_CO2, M_B, surfaceTension, wCO2_g_val,
129  wBA_g_val, prsr, bubble_radius
130  );
131  wCO2_Max =
132  wCO2Max
133  (
134  M_CO2, M_liq, partialPressure_CO2, henry_Coefficient
135  );
136 
137  if (blowingAgent == "n-pentane")
138  {
139  LMAX = LliqMax (tmp_val);
140  }
141  if (blowingAgent == "R-11")
142  {
143  xBL_val = xBL(tmp_val, dxdT);
144  LMAX = wBL_D(xBL_val, M_B, M_NCO, L0);
145  }
146  if (blowingAgent == "no")
147  {
148  LMAX = 0.0;
149  }
150 
151  if (bubbleGrowthSurrogateSwitch)
152  {
153  // Required variable for bubble growth surrogate models
154  double p_1 =
156  (
157  M_B, M_CO2, surfaceTension, wBA_g_val, wCO2_g_val,
158  prsr, bubble_radius
159  );
160  double p_2 = partialPressure_CO2;
161  double c_1 = wBA_l[celli]*rhoPolyS*1000.0/M_B;
162  double c_2 = wCO2_l[celli]*rhoPolyS*1000.0/M_CO2;
163  if (bubbleGrowthMode == "twoNodes")
164  {
165  // Info<< "p_1: " << p_1 << endl;
166  // Info<< "p_2: " << p_2 << endl;
167  // Info<< "c_1: " << c_1 << endl;
168  // Info<< "c_2: " << c_2 << endl;
169  // Info<< "bubble_radius: " << bubble_radius << endl;
170  // Info<< "T: " << TS[celli] << endl;
171  for (int i = 0; i < nodes; i++)
172  {
173  nodeRadii[i] = nodeRadius(vi[i]);
174  // set input vector
175  modena_inputs_set(inputs_bblgr1, Tbblgr1pos, TS[celli]);
176  modena_inputs_set(inputs_bblgr1, Rbblgr1pos, nodeRadii[i]);
177  // modena_inputs_set(inputs_bblgr1, KH1bblgr1pos, KH1);
178  modena_inputs_set(inputs_bblgr1, c_1bblgr1pos, c_1);
179  modena_inputs_set(inputs_bblgr1, p_1bblgr1pos, p_1);
180  // call the bblgr1 model
181  ret_bblgr1[i] =
182  modena_model_call (bblgr1, inputs_bblgr1, outputs_bblgr1);
183  if (ret_bblgr1[i] != 0)
184  {
185  modena_inputs_destroy (inputs_bblgr1);
186  modena_outputs_destroy (outputs_bblgr1);
187  modena_model_destroy (bblgr1);
188  exit(ret_bblgr1[i]);
189  }
190  radiusGrowthBA[i] =
191  modena_outputs_get(outputs_bblgr1, 0);
192 
193  // set input vector
194  modena_inputs_set(inputs_bblgr2, Tbblgr2pos, TS[celli]);
195  modena_inputs_set(inputs_bblgr2, Rbblgr2pos, nodeRadii[i]);
196  // modena_inputs_set(inputs_bblgr2, KH2bblgr2pos, KH2);
197  modena_inputs_set(inputs_bblgr2, c_2bblgr2pos, c_2);
198  modena_inputs_set(inputs_bblgr2, p_2bblgr2pos, p_2);
199  // call the bblgr2 model
200  ret_bblgr2[i] =
201  modena_model_call (bblgr2, inputs_bblgr2, outputs_bblgr2);
202  if (ret_bblgr2[i] != 0)
203  {
204  modena_inputs_destroy (inputs_bblgr2);
205  modena_outputs_destroy (outputs_bblgr2);
206  modena_model_destroy (bblgr2);
207  exit(ret_bblgr2[i]);
208  }
209  radiusGrowthCO2[i] =
210  modena_outputs_get(outputs_bblgr2, 0);
211 
212  volumeGrowthBA[i] =
213  (radiusGrowthBA[i]*RR*TS[celli])/(p_1);
214  volumeGrowthCO2[i] =
215  (radiusGrowthCO2[i]*RR*TS[celli])/(p_2);
216  // Info<< "radiusGrowthBA[" << i << "] " << radiusGrowthBA[i] << endl;
217  // Info<< "radiusGrowthCO2[" << i << "] " << radiusGrowthCO2[i] << endl;
218  // Info<< "volumeGrowthBA[" << i << "] " << volumeGrowthBA[i] << endl;
219  // Info<< "volumeGrowthCO2[" << i << "] " << volumeGrowthCO2[i] << endl;
220  if
221  (
222  volumeGrowthBA[i] < 0.0
223  || radiusGrowthBA[i] < 0.0
224  || L0 < 1.0e-8
225  || wBA_l[celli] < 1.0e-8
226  )
227  {
228  volumeGrowthBA[i] = 0.0;
229  }
230  if
231  (
232  volumeGrowthCO2[i] < 0.0
233  || radiusGrowthCO2[i] < 0.0
234  || CW_0 < 1.0e-8
235  || wCO2_l[celli] < 1.0e-8
236  )
237  {
238  volumeGrowthCO2[i] = 0.0;
239  }
240  }
241  }
242  else if (bubbleGrowthMode == "meanRadius")
243  {
244  // Info<< "p_1: " << p_1 << endl;
245  // Info<< "p_2: " << p_2 << endl;
246  // Info<< "c_1: " << c_1 << endl;
247  // Info<< "c_2: " << c_2 << endl;
248  // Info<< "bubble_radius: " << bubble_radius << endl;
249  // Info<< "T: " << TS[celli] << endl;
250 
251  // set input vector
252  modena_inputs_set(inputs_bblgr1, Tbblgr1pos, TS[celli]);
253  modena_inputs_set(inputs_bblgr1, Rbblgr1pos, bubble_radius);
254  modena_inputs_set(inputs_bblgr1, c_1bblgr1pos, c_1);
255  modena_inputs_set(inputs_bblgr1, p_1bblgr1pos, p_1);
256  // set input vector
257  modena_inputs_set(inputs_bblgr2, Tbblgr2pos, TS[celli]);
258  modena_inputs_set(inputs_bblgr2, Rbblgr2pos, bubble_radius);
259  modena_inputs_set(inputs_bblgr2, c_2bblgr2pos, c_2);
260  modena_inputs_set(inputs_bblgr2, p_2bblgr2pos, p_2);
261  // call the bblgr1 model
262  int ret_bblgr_1 = modena_model_call (bblgr1, inputs_bblgr1, outputs_bblgr1);
263  if (ret_bblgr_1 != 0)
264  {
265  modena_inputs_destroy (inputs_bblgr1);
266  modena_outputs_destroy (outputs_bblgr1);
267  modena_model_destroy (bblgr1);
268  exit(ret_bblgr_1);
269  }
270  // call the bblgr2 model
271  int ret_bblgr_2 = modena_model_call (bblgr2, inputs_bblgr2, outputs_bblgr2);
272  if (ret_bblgr_2 != 0)
273  {
274  modena_inputs_destroy (inputs_bblgr2);
275  modena_outputs_destroy (outputs_bblgr2);
276  modena_model_destroy (bblgr2);
277  exit(ret_bblgr_2);
278  }
279  double G1, G2;
280  for (int i = 0; i < nNodes; i++)
281  {
282  volumeGrowthBA[i] = 0.0;
283  volumeGrowthCO2[i] = 0.0;
284  }
285  G1 = modena_outputs_get(outputs_bblgr1, 0);
286  G2 = modena_outputs_get(outputs_bblgr2, 0);
287  // Info<< "G1: " << G1 << endl;
288  // Info<< "G2: " << G2 << endl;
289 
290  volumeGrowthBA[0] = (G1*RR*TS[celli])/(p_1);
291  volumeGrowthCO2[0] = (G2*RR*TS[celli])/(p_2);
292  // Info<< "volumeGrowthBA[0]" << volumeGrowthBA[0] << endl;
293  // Info<< "volumeGrowthCO2[0]" << volumeGrowthCO2[0] << endl;
294  if
295  (
296  volumeGrowthBA[0] < 0.0
297  || G1 < 0.0
298  || L0 < 1.0e-8
299  || wBA_l[celli] < 1.0e-8
300  )
301  {
302  volumeGrowthBA[0] = 0.0;
303  }
304  if
305  (
306  volumeGrowthCO2[0] < 0.0
307  || G2 < 0.0
308  || CW_0 < 1.0e-8
309  || wCO2_l[celli] < 1.0e-8
310  )
311  {
312  volumeGrowthCO2[0] = 0.0;
313  }
314  volumeGrowthBA[1] = volumeGrowthBA[0];
315  volumeGrowthCO2[1] = volumeGrowthCO2[0];
316  }
317  else
318  {
319  FatalErrorIn(args.executable())
320  << "Invalid bubble growth mode " << bubbleGrowthMode
321  << "\nValid modes are: "
322  << "\ntwoNodes,"
323  << "\nmeanRadius"
324  << exit(FatalError);
325  }
326  }
327  else
328  {
329  for (int i = 0; i < nodes; i++)
330  {
331  volumeGrowthBA[i] = growthRateBA(L_l_val,LMAX); // growth rate constant due to blowing agent
332  volumeGrowthCO2[i]= growthRateCO2(CO2_l_val, wCO2_Max); // growth rate constant due to CO2
333  }
334  }
335 
337  (
338  sgBA, sgCO2, we, vi, nodes, mOrder, CO2_l_val, L_l_val,
339  tmp_val, wCO2_Max, cc1_val, LMAX,
340  volumeGrowthBA, volumeGrowthCO2
341  );
342 
343  g1_CO2[celli] = sgCO2[1];
344  g1_BA[celli] = sgBA[1];
345  cc1[celli] = cc1_val;
346 
347  if (XOH[celli] > XOH_Gel)
348  {
349  for (int j = 0; j < 2*nodes; j++)
350  {
351  sc[j] = 0.0;
352  }
353  }
354  else
355  {
356  coalescenceSource(sc, we, vi, nodes, mOrder);
357  }
358 
359  for (int j = 0; j < 2*nodes; j++)
360  {
361  sources[j] = sgBA[j] + sgCO2[j] + 0.5*sc[j];
362  }
363  mZeroSource[celli] = sources[0];
364  mOneSource[celli] = sources[1];
365  mTwoSource[celli] = sources[2];
366  mThreeSource[celli] = sources[3];
367  }
368  else
369  {
370  mZeroSource[celli] = ROOTVSMALL;
371  mOneSource[celli] = ROOTVSMALL;
372  mTwoSource[celli] = ROOTVSMALL;
373  mThreeSource[celli] = ROOTVSMALL;
374  }
375  } //end of forAll
376  }
377 
378  if
379  (
380  PBEMethod == "QMOM"
381  && nNodes == 3
382  )
383  {
384  static int nodes = nNodes;
385  int mOrder[2*nodes];
386  for (int i = 0; i < 2*nNodes; i++){mOrder[i] = i;}
387  double mom[2*nodes], we[nodes], vi[nodes], sgBA[2*nodes],
388  sgCO2[2*nodes], sc[2*nodes];
389  double sources[2*nodes];
390 
391  double CO2_l_val, L_l_val, tmp_val, alphafoam;
392  double LMAX, xBL_val, rhoPolyS;
393  double cc1_val = 0.0;
394 
395  double radiusGrowthBA[nodes],
396  radiusGrowthCO2[nodes],
397  volumeGrowthBA[nodes],
398  volumeGrowthCO2[nodes],
399  nodeRadii[nodes];
400  int ret_bblgr1[nodes], ret_bblgr2[nodes];
401  // Calling the model for density reaction mixture
402  size_t Td_pos =
403  modena_model_inputs_argPos(density_reaction_mixturemodel, "T");
404  size_t XOHd_pos =
405  modena_model_inputs_argPos(density_reaction_mixturemodel, "XOH");
406  modena_model_argPos_check(density_reaction_mixturemodel);
407 
408  size_t Tbblgr1pos = modena_model_inputs_argPos(bblgr1, "T");
409  size_t Rbblgr1pos = modena_model_inputs_argPos(bblgr1, "R");
410  // size_t KH1bblgr1pos = modena_model_inputs_argPos(bblgr1, "kH");
411  size_t c_1bblgr1pos = modena_model_inputs_argPos(bblgr1, "c");
412  size_t p_1bblgr1pos = modena_model_inputs_argPos(bblgr1, "p");
414  size_t Tbblgr2pos = modena_model_inputs_argPos(bblgr2, "T");
415  size_t Rbblgr2pos = modena_model_inputs_argPos(bblgr2, "R");
416  // size_t KH2bblgr2pos = modena_model_inputs_argPos(bblgr2, "kH");
417  size_t c_2bblgr2pos = modena_model_inputs_argPos(bblgr2, "c");
418  size_t p_2bblgr2pos = modena_model_inputs_argPos(bblgr2, "p");
420 
421  forAll(mesh.C(), celli)
422  {
423  CO2_l_val = wCO2_l[celli];
424  L_l_val = wBA_l[celli];
425  tmp_val = TS[celli];
426  alphafoam = alpha2[celli];
427 
428  // Required variables for wCO2_Max
429  double wCO2_g_val, wBA_g_val,prsr;
430  wCO2_g_val= wCO2_g[celli];
431  wBA_g_val = wBA_g[celli];
432  prsr = p[celli];
433 
434  if (liquidMixtureDensitySurrogate)
435  {
436  // set input vector
437  modena_inputs_set(inputs_den, Td_pos, TS[celli]);
438  modena_inputs_set(inputs_den, XOHd_pos, XOH[celli]);
439  int ret_den =
441  (
442  density_reaction_mixturemodel,
443  inputs_den,
444  outputs_den
445  );
446  if (ret_den != 0)
447  {
448  modena_inputs_destroy (inputs_den);
449  modena_outputs_destroy (outputs_den);
450  modena_model_destroy (density_reaction_mixturemodel);
451  exit(ret_den);
452  }
453  rhoPolyS = modena_outputs_get(outputs_den, 0);
454  }
455  else
456  {
457  rhoPolyS = rhoPoly;
458  }
459 
460  if (alphafoam > 0.5)
461  {
462  mom[0] = mZero[celli];
463  mom[1] = mOne[celli];
464  mom[2] = mTwo[celli];
465  mom[3] = mThree[celli];
466  mom[4] = mFour[celli];
467  mom[5] = mFive[celli];
468 
469  // calling the functions, PDA, growthSource, coalescenceSource
470  PDA(we, vi, mom, nodes);
471  weight0[celli] = we[0];
472  weight1[celli] = we[1];
473  weight2[celli] = we[2];
474  node0[celli] = vi[0];
475  node1[celli] = vi[1];
476  node2[celli] = vi[2];
477 
478  // Required variables for wCO2_Max
479  double henry_Coefficient, bubble_radius,
480  partialPressure_CO2, wCO2_Max;
481  henry_Coefficient = henryCoefficient(tmp_val);
482  bubble_radius = bubbleRadius(mom[0], mom[1]);
483  partialPressure_CO2 =
485  (
486  M_CO2, M_B, surfaceTension, wCO2_g_val, wBA_g_val,
487  prsr, bubble_radius
488  );
489  wCO2_Max =
490  wCO2Max
491  (
492  M_CO2, M_liq, partialPressure_CO2, henry_Coefficient
493  );
494 
495  if (blowingAgent == "n-pentane")
496  {
497  LMAX = LliqMax (tmp_val);
498  }
499  if (blowingAgent == "R-11")
500  {
501  xBL_val = xBL(tmp_val, dxdT);
502  LMAX = wBL_D(xBL_val, M_B, M_NCO, L0);
503  }
504  if (blowingAgent == "no")
505  {
506  LMAX = 0.0;
507  }
508 
509  if (bubbleGrowthSurrogateSwitch)
510  {
511  // Required variable for bubble growth surrogate models
512  double p_1 =
514  (
515  M_B, M_CO2, surfaceTension, wBA_g_val,
516  wCO2_g_val, prsr, bubble_radius
517  );
518  double p_2 = partialPressure_CO2;
519  double c_1 = wBA_l[celli]*rhoPolyS*1000.0/M_B;
520  double c_2 = wCO2_l[celli]*rhoPolyS*1000.0/M_CO2;
521 
522  for (int i = 0; i < nodes; i++)
523  {
524  nodeRadii[i] = nodeRadius(vi[i]);
525  // set input vector
526  modena_inputs_set(inputs_bblgr1, Tbblgr1pos, TS[celli]);
527  modena_inputs_set(inputs_bblgr1, Rbblgr1pos, nodeRadii[i]);
528  // modena_inputs_set(inputs_bblgr1, KH1bblgr1pos, KH1);
529  modena_inputs_set(inputs_bblgr1, c_1bblgr1pos, c_1);
530  modena_inputs_set(inputs_bblgr1, p_1bblgr1pos, p_1);
531  // call the bblgr1 model
532  ret_bblgr1[i] =
533  modena_model_call (bblgr1, inputs_bblgr1, outputs_bblgr1);
534  if (ret_bblgr1[i] != 0)
535  {
536  modena_inputs_destroy (inputs_bblgr1);
537  modena_outputs_destroy (outputs_bblgr1);
538  modena_model_destroy (bblgr1);
539  exit(ret_bblgr1[i]);
540  }
541  radiusGrowthBA[i] =
542  modena_outputs_get(outputs_bblgr1, 0);
543 
544  // set input vector
545  modena_inputs_set(inputs_bblgr2, Tbblgr2pos, TS[celli]);
546  modena_inputs_set(inputs_bblgr2, Rbblgr2pos, nodeRadii[i]);
547  // modena_inputs_set(inputs_bblgr2, KH2bblgr2pos, KH2);
548  modena_inputs_set(inputs_bblgr2, c_2bblgr2pos, c_2);
549  modena_inputs_set(inputs_bblgr2, p_2bblgr2pos, p_2);
550  // call the bblgr2 model
551  ret_bblgr2[i] =
552  modena_model_call (bblgr2, inputs_bblgr2, outputs_bblgr2);
553  if (ret_bblgr2[i] != 0)
554  {
555  modena_inputs_destroy (inputs_bblgr2);
556  modena_outputs_destroy (outputs_bblgr2);
557  modena_model_destroy (bblgr2);
558  exit(ret_bblgr2[i]);
559  }
560  radiusGrowthCO2[i] =
561  modena_outputs_get (outputs_bblgr2, 0);
562 
563  volumeGrowthBA[i] =
564  (radiusGrowthBA[i]*RR*TS[celli])/(p_1);
565  volumeGrowthCO2[i] =
566  (radiusGrowthCO2[i]*RR*TS[celli])/(p_2);
567  if
568  (
569  volumeGrowthBA[i] < 0.0
570  || radiusGrowthBA[i] < 0.0
571  || L0 < 1.0e-8
572  || wBA_l[celli] < 1.0e-8
573  )
574  {
575  volumeGrowthBA[i] = 0.0;
576  }
577  if
578  (
579  volumeGrowthCO2[i] < 0.0
580  || radiusGrowthCO2[i] < 0.0
581  || CW_0 < 1.0e-8
582  || wCO2_l[celli] < 1.0e-8
583  )
584  {
585  volumeGrowthCO2[i] = 0.0;
586  }
587  }
588  }
589  else
590  {
591  for (int i = 0; i < nodes; i++)
592  {
593  volumeGrowthBA[i] = growthRateBA(L_l_val,LMAX); // growth rate constant due to blowing agent
594  volumeGrowthCO2[i]= growthRateCO2(CO2_l_val, wCO2_Max); // growth rate constant due to CO2
595  }
596  }
598  (
599  sgBA, sgCO2, we, vi, nodes, mOrder, CO2_l_val, L_l_val,
600  tmp_val, wCO2_Max, cc1_val, LMAX,
601  volumeGrowthBA, volumeGrowthCO2
602  );
603 
604  g1_CO2[celli] = sgCO2[1];
605  g1_BA[celli] = sgBA[1];
606  cc1[celli] = cc1_val;
607 
608  if (XOH[celli] > XOH_Gel)
609  {
610  for (int j = 0; j < 2*nodes; j++)
611  {
612  sc[j] = 0.0;
613  }
614  }
615  else
616  {
617  coalescenceSource(sc, we, vi, nodes, mOrder);
618  }
619 
620  for (int j = 0; j < 2*nodes; j++)
621  {
622  sources[j] = sgBA[j] + sgCO2[j] + 0.5*sc[j];
623  }
624  mZeroSource[celli] = sources[0];
625  mOneSource[celli] = sources[1];
626  mTwoSource[celli] = sources[2];
627  mThreeSource[celli] = sources[3];
628  mFourSource[celli] = sources[4];
629  mFiveSource[celli] = sources[5];
630  }
631  else
632  {
633  mZeroSource[celli] = ROOTVSMALL;
634  mOneSource[celli] = ROOTVSMALL;
635  mTwoSource[celli] = ROOTVSMALL;
636  mThreeSource[celli] = ROOTVSMALL;
637  mFourSource[celli] = ROOTVSMALL;
638  mFiveSource[celli] = ROOTVSMALL;
639  }
640  } //end of forAll
641  }
642 }
643 ///@endcond
size_t modena_model_inputs_argPos(const modena_model_t *self, const char *name)
Function determining position of an argument in the input vector.
Definition: model.c:366
double henryCoefficient(double &T)
Henry coefficient for CO2.
void growthSource(double *, double *, double *, double *, int &, int *, double &, double &, double &, double *, double *)
Definition: growth.h:22
void coalescenceSource(double *, double *, double *, int &, int *, double &)
Definition: coalescence.h:17
double partialPressureCO2(const state_type &y)
partial pressure of CO2
double surfaceTension
required for the computation of partial pressure
double M_CO2
Molecular mass of carbon dioxide, kg/kmol.
double M_NCO
Molecular weight of NCO, kg/kmol.
int modena_model_call(modena_model_t *self, modena_inputs_t *inputs, modena_outputs_t *outputs)
Function calling the surrogate model and checking for errors.
Definition: model.c:553
double RR
ideal gas constant, J/mol K
double bubbleRadius(const double m0, const double m1)
radius of bubbles based on the moments
Definition: bubbleRadius.h:28
double LliqMax(double &)
maximum allowable amount of liquid blowing agent (n-pentane) in liquid
double wBL_D(double &xBL, double &M_B, double &M_NCO, double &L0)
weight fraction of maximum allowable blowing agent (R-11) in liquid
double rhoPoly
density of the liquid polymer, kg/m3
double growthRateBA(double &, double &)
growth rate due to evaporation of blowing agents
Definition: growthRate.H:41
void modena_model_destroy(modena_model_t *self)
Function deallocating the memory allocated for the surrogate model.
Definition: model.c:665
double nodeRadius(const double &v)
radius of bubbles at each node
Definition: bubbleRadius.h:38
double M_B
Molecular mass of blowing agent, kg/kmol.
double xBL(double &T, double &dxdT)
mole fraction of blowing agent (R-11) in liquid polymer
double partialPressureBA(const state_type &y)
partial pressure of the physical blowing agent
double L0
Initial weight fraction of blowing agent in the liquid, -.
double wCO2Max(double &M_CO2, double &M_liq, double &pCO2, double &henryCoeff)
dissolved amount of CO2 in liquid
void PDA(double *, double *, double *, int &)
Definition: pda.h:13
void modena_model_argPos_check(const modena_model_t *self)
Function checking that the user has queried all input positions.
Definition: model.c:408
double growthRateCO2(double &, double &)
growth rate due to production of CO2
Definition: growthRate.H:28