TrioCFD 1.9.9_beta
TrioCFD documentation
Loading...
Searching...
No Matches
Op_Diff_VEF_Face_Stab.cpp
1/****************************************************************************
2* Copyright (c) 2024, CEA
3* All rights reserved.
4*
5* Redistribution and use in source and binary forms, with or without modification, are permitted provided that the following conditions are met:
6* 1. Redistributions of source code must retain the above copyright notice, this list of conditions and the following disclaimer.
7* 2. Redistributions in binary form must reproduce the above copyright notice, this list of conditions and the following disclaimer in the documentation and/or other materials provided with the distribution.
8* 3. Neither the name of the copyright holder nor the names of its contributors may be used to endorse or promote products derived from this software without specific prior written permission.
9*
10* THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE DISCLAIMED.
11* IN NO EVENT SHALL THE COPYRIGHT HOLDER OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS;
12* OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
13*
14*****************************************************************************/
15
16#include <Op_Diff_VEF_Face_Stab.h>
17#include <Champ_P1NC.h>
18#include <Champ_Q1NC.h>
19
20#include <Periodique.h>
21#include <Symetrie.h>
22#include <Neumann_paroi.h>
23#include <Echange_externe_impose.h>
24#include <Neumann_sortie_libre.h>
25#include <Dirichlet.h>
26#include <Dirichlet_homogene.h>
27#include <Neumann_homogene.h>
28
29
30#include <Milieu_base.h>
31
32#include <TRUSTTrav.h>
33#include <Probleme_base.h>
34#include <Navier_Stokes_std.h>
35#include <Porosites_champ.h>
36
37#include <SFichier.h>
38#include <ArrOfBit.h>
39#include <Schema_Temps_base.h>
40
41Implemente_instanciable_sans_constructeur(Op_Diff_VEF_Face_Stab,"Op_Diff_VEFSTAB_P1NC",Op_Diff_VEF_Face);
42
43double minimum(double a,double b,double c)
44{
45 if (a<=b)
46 if (a<=c) return a;
47 else return c;
48 else if (b<=c) return b;
49 else return c;
50}
51
52double maximum(double a,double b,double c)
53{
54 if (a>=b)
55 if (a>=c) return a;
56 else return c;
57 else if (b>=c) return b;
58 else return c;
59}
60
61double minimum(double a,double b)
62{
63 if (a<=b) return a;
64 else return b;
65}
66
67double maximum(double a,double b)
68{
69 if (a>=b) return a;
70 else return b;
71}
72
77
79{
80 return s << que_suis_je() ;
81}
82
84{
85 Motcle motlu,accouverte="{",accfermee="}";
86 Motcles les_mots(3);
87 {
88 les_mots[0]="standard";
89 les_mots[1]="info";
90 les_mots[2]="new_jacobian";
91 }
92
93 s>>motlu;
94 if (motlu!=accouverte)
95 {
96 Cerr<<"Error in Op_Diff_VEF_Face_Stab::readOn()"<<finl;
97 Cerr<<"Option keywords must be preceded by an open brace"<<finl;
99 }
100
101 s>>motlu;
102 while (motlu!=accfermee)
103 {
104 int rang=les_mots.search(motlu);
105 switch(rang)
106 {
107 case 0 :
108 s>>standard_;
109 break;
110 case 1 :
111 s>>info_;
112 break;
113 case 2 :
114 s>>new_jacobian_;
115 break;
116 default :
117 Cerr<<"Error in Op_Diff_VEF_Face_Stab::readOn()"<<finl;
118 Cerr<<"Word "<<motlu<<" is unknown"<<finl;
119 Cerr<<"Known keywords are : "<<les_mots<<finl;
121 break;
122 }
123 s>>motlu;
124 }
125 return s;
126}
127
128
129DoubleTab& Op_Diff_VEF_Face_Stab::ajouter(const DoubleTab& inconnue_org, DoubleTab& resu) const
130{
132 const Domaine_VEF& domaine_VEF = le_dom_vef.valeur();
133
134 const int nb_elem_tot=domaine_VEF.nb_elem_tot();
135 const int nb_faces_elem=domaine_VEF.domaine().nb_faces_elem();
136
137 // either div(phi nu grad inco) or div(nu grad phi inco)
138 // depending on whether phi_psi or psi is being diffused
139 DoubleTab nu;
140 DoubleTab tab_inconnue;
141
142 int marq=phi_psi_diffuse(equation());
143
144 const DoubleVect& porosite_face = equation().milieu().porosite_face();
145 const DoubleVect& porosite_elem = equation().milieu().porosite_elem();
146
147 modif_par_porosite_si_flag(nu_,nu,!marq,porosite_elem);
148 const DoubleTab& inconnue=modif_par_porosite_si_flag(inconnue_org,tab_inconnue,marq,porosite_face);
149 //
150 //
151 //
152 DoubleTab resu2(resu);
153 resu2=0.;
154
155 DoubleTab Aij(nb_elem_tot,nb_faces_elem,nb_faces_elem);
156 calculer_coefficients(nu,Aij);
157
158 ajouter_operateur_centre(Aij,inconnue,resu2);
159 if (!standard_)
160 {
161 ajouter_diffusion(Aij,inconnue,resu2);
162 ajouter_antidiffusion(Aij,inconnue,resu2);
163 }
164 modifie_pour_Cl(inconnue,resu2);
165
166 modifier_flux(*this);
167
168 resu-=resu2;//-= because the Laplacian is placed as a source term in the equation
169 return resu;
170}
171
172void Op_Diff_VEF_Face_Stab::modifie_pour_Cl(const DoubleTab& inconnue, DoubleTab& resu) const
173{
174 const Domaine_VEF& domaine_VEF = le_dom_vef.valeur();
175 const Domaine_Cl_VEF& domaine_Cl_VEF = la_zcl_vef.valeur();
176 const Conds_lim& les_cl = domaine_Cl_VEF.les_conditions_limites();
177
178 const DoubleVect& inconnueVect = inconnue;
179 DoubleVect& resuVect = resu;
180
181 const int nb_bords=les_cl.size();
182
183 int nb_comp=1;
184 if(resu.nb_dim()==2) nb_comp=resu.dimension(1);
185
186 DoubleTab& tab_flux_bords = flux_bords_;
187 tab_flux_bords.resize(domaine_VEF.nb_faces_bord(),nb_comp);
188 tab_flux_bords=0.;
189
190 int num1=0, num2=0;
191 int n_bord=0;
192 int face=0;
193 int face_associee=0;
194 int ind_face=0;
195 int dim=0;
196
197 double surface=0.;
198 double flux=0.;
199
200 for (n_bord=0; n_bord<nb_bords; n_bord++)
201 {
202 const Cond_lim& la_cl = domaine_Cl_VEF.les_conditions_limites(n_bord);
203 const Front_VF& le_bord = ref_cast(Front_VF,la_cl->frontiere_dis());
204
205 num1 = 0;
206 num2 = le_bord.nb_faces();
207
208 if (sub_type(Periodique,la_cl.valeur()))
209 {
210 const Periodique& la_cl_perio = ref_cast(Periodique,la_cl.valeur());
211
212 for (ind_face=num1; ind_face<num2; ind_face++)
213 {
214 face=le_bord.num_face(ind_face);
215 face_associee=le_bord.num_face(la_cl_perio.face_associee(ind_face));
216
217 if (face<face_associee)
218 for (dim=0; dim<nb_comp; dim++)
219 {
220 resuVect[face*nb_comp+dim]+=resuVect[face_associee*nb_comp+dim];
221 resuVect[face_associee*nb_comp+dim]=resuVect[face*nb_comp+dim];
222 }
223 }
224 }
225 if (sub_type(Neumann_paroi,la_cl.valeur()))
226 {
227 const Neumann_paroi& la_cl_paroi = ref_cast(Neumann_paroi,la_cl.valeur());
228
229 for (ind_face=num1; ind_face<num2; ind_face++)
230 {
231 face=le_bord.num_face(ind_face);
232 surface=domaine_VEF.surface(face);
233
234 for (dim=0; dim<nb_comp; dim++)
235 {
236 flux=la_cl_paroi.flux_impose(ind_face,dim)*surface;
237 resuVect[face*nb_comp+dim]-=flux;
238 tab_flux_bords(face,dim)=flux;
239 }
240 }
241 }
242 if (sub_type(Echange_externe_impose,la_cl.valeur()))
243 {
244 const Echange_externe_impose& la_cl_paroi = ref_cast(Echange_externe_impose,la_cl.valeur());
245
246 for (ind_face=num1; ind_face<num2; ind_face++)
247 {
248 face=le_bord.num_face(ind_face);
249 surface=domaine_VEF.surface(face);
250
251 for (dim=0; dim<nb_comp; dim++)
252 {
253 flux=la_cl_paroi.h_imp(ind_face,dim)*(la_cl_paroi.T_ext(ind_face,dim)-inconnueVect[face*nb_comp+dim])*surface;
254 resuVect[face*nb_comp+dim]-=flux;
255 tab_flux_bords(face,dim)=flux;
256 }
257 }
258 }
259 if (sub_type(Neumann_homogene,la_cl.valeur())
260 || sub_type(Symetrie,la_cl.valeur())
261 || sub_type(Neumann_sortie_libre,la_cl.valeur()))
262 {
263 for (ind_face=num1; ind_face<num2; ind_face++)
264 {
265 face=le_bord.num_face(ind_face);
266
267 for (dim=0; dim<nb_comp; dim++)
268 tab_flux_bords(face,dim)=0.;
269 }
270 }
271
272 }
273}
274
275void Op_Diff_VEF_Face_Stab::ajouter_operateur_centre(const DoubleTab& Aij, const DoubleTab& inconnueTab, DoubleTab& resuTab) const
276{
277 const Domaine_VEF& domaine_VEF = le_dom_vef.valeur();
278
279 const int nb_elem_tot=domaine_VEF.nb_elem_tot();
280 const int nb_faces_elem=domaine_VEF.domaine().nb_faces_elem();
281
282 int nb_comp=1;
283 if(resuTab.nb_dim()==2) nb_comp=resuTab.dimension(1);
284
285 const DoubleVect& inconnue = inconnueTab;
286 DoubleVect& resu = resuTab;
287
288 int elem=0;
289 int facei_loc=0,facei=0;
290 int facej_loc=0,facej=0;
291 int dim=0;
292
293 double inc_i=0.;
294 double inc_j=0.;
295 double delta_ij=0.;
296
297 const IntTab& elem_faces=domaine_VEF.elem_faces();
298
299 for (elem=0; elem<nb_elem_tot; elem++)
300 for (facei_loc=0; facei_loc<nb_faces_elem; facei_loc++)
301 {
302 facei=elem_faces(elem,facei_loc);
303
304 for (facej_loc=facei_loc+1; facej_loc<nb_faces_elem; facej_loc++)
305 {
306 facej=elem_faces(elem,facej_loc);
307
308 const double aij=Aij(elem,facei_loc,facej_loc);
309
310 for (dim=0; dim<nb_comp; dim++)
311 {
312 inc_i=inconnue[facei*nb_comp+dim];
313 inc_j=inconnue[facej*nb_comp+dim];
314 delta_ij=aij*(inc_j-inc_i);
315
316 resu[facei*nb_comp+dim]+=delta_ij;
317 resu[facej*nb_comp+dim]-=delta_ij;
318 }
319 }
320 }
321}
322
323void Op_Diff_VEF_Face_Stab::ajouter_diffusion(const DoubleTab& Aij, const DoubleTab& inconnueTab, DoubleTab& resuTab) const
324{
325 const Domaine_VEF& domaine_VEF = le_dom_vef.valeur();
326
327 const int nb_elem_tot=domaine_VEF.nb_elem_tot();
328 const int nb_faces_elem=domaine_VEF.domaine().nb_faces_elem();
329
330 int nb_comp=1;
331 if(resuTab.nb_dim()==2) nb_comp=resuTab.dimension(1);
332
333 const DoubleVect& inconnue= inconnueTab;
334 DoubleVect& resu= resuTab;
335
336 int elem=0;
337 int facei_loc=0,facei=0;
338 int facej_loc=0,facej=0;
339 int dim=0;
340
341 double inc_i=0.;
342 double inc_j=0.;
343 double delta_ij=0.;
344 double dij=0.;
345
346 const IntTab& elem_faces=domaine_VEF.elem_faces();
347
348 for (elem=0; elem<nb_elem_tot; elem++)
349 for (facei_loc=0; facei_loc<nb_faces_elem; facei_loc++)
350 {
351 facei=elem_faces(elem,facei_loc);
352
353 for (facej_loc=facei_loc+1; facej_loc<nb_faces_elem; facej_loc++)
354 {
355 facej=elem_faces(elem,facej_loc);
356
357 const double aij=Aij(elem,facei_loc,facej_loc);
358
359 if (aij>0.)
360 {
361 dij=-aij;
362
363 for (dim=0; dim<nb_comp; dim++)
364 {
365 inc_i=inconnue[facei*nb_comp+dim];
366 inc_j=inconnue[facej*nb_comp+dim];
367 delta_ij=dij*(inc_j-inc_i);
368
369 resu[facei*nb_comp+dim]+=delta_ij;
370 resu[facej*nb_comp+dim]-=delta_ij;
371 }
372 }
373 }
374 }
375}
376
377void Op_Diff_VEF_Face_Stab::ajouter_antidiffusion(const DoubleTab& Aij, const DoubleTab& inconnueTab, DoubleTab& resuTab) const
378{
379 const Domaine_VEF& domaine_VEF = le_dom_vef.valeur();
380
381 const int nb_elem_tot=domaine_VEF.nb_elem_tot();
382 const int nb_faces_tot=domaine_VEF.nb_faces_tot();
383 const int nb_faces_elem=domaine_VEF.domaine().nb_faces_elem();
384
385 int nb_comp=1;
386 if(resuTab.nb_dim()==2) nb_comp=resuTab.dimension(1);
387
388 const DoubleVect& inconnue = inconnueTab;
389 DoubleVect& resu = resuTab;
390
391 const DoubleTab& xv=domaine_VEF.xv();
392 DoubleTab rij(Objet_U::dimension);
393 DoubleTab rji(rij);
394
395 int elem=0;
396 int facei_loc=0,facei=0;
397 int facej_loc=0,facej=0;
398 int dim=0;
399 int dim2=0;
400
401 double inc_i=0.;
402 double inc_j=0.;
403 double delta_ij=0.;
404 double delta_imax=0.,delta_imin=0.;
405 double delta_jmax=0.,delta_jmin=0.;
406 double sij=0.;
407 double sij_max=DMINFLOAT;
408 double sij_min=DMAXFLOAT;
409 double muij=0.,muji=0.;
410
411 bool ok_facei=false;
412 bool ok_facej=false;
413
414
415
416 const IntTab& elem_faces=domaine_VEF.elem_faces();
417
418 DoubleTab Minima(nb_faces_tot);
419 DoubleTab Maxima(nb_faces_tot);
420
421
422 for (dim=0; dim<nb_comp; dim++)
423 {
424 calculer_min(inconnueTab,dim,Minima);
425 calculer_max(inconnueTab,dim,Maxima);
426
427 for (elem=0; elem<nb_elem_tot; elem++)
428 for (facei_loc=0; facei_loc<nb_faces_elem; facei_loc++)
429 {
430 facei=elem_faces(elem,facei_loc);
431 ok_facei=is_dirichlet_faces_(facei);
432 inc_i=inconnue[facei*nb_comp+dim];
433
434 delta_imin=Minima(facei)-inc_i;
435 delta_imax=Maxima(facei)-inc_i;
436 assert(delta_imax>=0.);
437 assert(delta_imin<=0.);
438
439 for (facej_loc=facei_loc+1; facej_loc<nb_faces_elem; facej_loc++)
440 {
441 const double aij=Aij(elem,facei_loc,facej_loc);
442
443 if (aij>0.)
444 {
445 facej=elem_faces(elem,facej_loc);
446 ok_facej=is_dirichlet_faces_(facej);
447 inc_j=inconnue[facej*nb_comp+dim];
448
449 for (dim2=0; dim2<Objet_U::dimension; dim2++)
450 rij(dim2)=xv(facej,dim2)-xv(facei,dim2);
451
452 rji=rij;
453 rji*=-1.;
454
455 delta_ij=inc_i-inc_j;
456 delta_jmin=Minima(facej)-inc_j;
457 delta_jmax=Maxima(facej)-inc_j;
458 assert(delta_jmax>=0.);
459 assert(delta_jmin<=0.);
460
461 muij=calculer_gradients(facei,rij);
462 muji=calculer_gradients(facej,rji);
463
464 sij=0.; //remains 0 if all faces are Dirichlet
465 if (delta_ij>0.)
466 {
467 muij*=delta_imax;
468 muji*=-delta_jmin;
469
470 if (!ok_facei && !ok_facej) //no Dirichlet face
471 sij=minimum(muij,delta_ij,muji);
472 if (!ok_facei && ok_facej) //facej is Dirichlet, facei is not
473 sij=minimum(muij,delta_ij);
474 if (ok_facei && !ok_facej) //facei is Dirichlet, facej is not
475 sij=minimum(delta_ij,muji);
476 }
477 else if (delta_ij<0.)
478 {
479 muij*=delta_imin;
480 muji*=-delta_jmax;
481
482 if (!ok_facei && !ok_facej) //no Dirichlet face
483 sij=maximum(muij,delta_ij,muji);
484 if (!ok_facei && ok_facej) //facej is Dirichlet, facei is not
485 sij=maximum(muij,delta_ij);
486 if (ok_facei && !ok_facej) //facei is Dirichlet, facej is not
487 sij=maximum(delta_ij,muji);
488 }
489
490 if (info_)
491 if (delta_ij!=0.)
492 {
493 double tmp=sij/delta_ij;
494 if (tmp>sij_max) sij_max=tmp;
495 if (tmp<sij_min) sij_min=tmp;
496 }
497
498 resu[facei*nb_comp+dim]-=aij*sij;
499 resu[facej*nb_comp+dim]+=aij*sij;
500
501 }
502 }
503 }
504 }
505
506 {
507 sij_max = Process::mp_max(sij_max);
508 sij_min = Process::mp_min(sij_min);
509 if (info_ && Process::je_suis_maitre())
510 {
511 SFichier mem_fichier("sij_memory.txt");
512 mem_fichier<<sij_max<<" "<<sij_min<<finl;
513 }
514 }
515}
516
517void Op_Diff_VEF_Face_Stab::calculer_coefficients(const DoubleTab& nu, DoubleTab& Aij) const
518{
519 const Domaine_VEF& domaine_VEF = le_dom_vef.valeur();
520
521 const int nb_elem_tot=domaine_VEF.nb_elem_tot();
522 const int nb_faces_elem=domaine_VEF.domaine().nb_faces_elem();
523
524 const IntTab& elem_faces=domaine_VEF.elem_faces();
525 const IntTab& face_voisins=domaine_VEF.face_voisins();
526
527 const DoubleTab& face_normales=domaine_VEF.face_normales();
528 const DoubleVect& volumes=domaine_VEF.volumes();
529
530 double volume=0.;
531 double signei=0.;
532 double signej=0.;
533 double psc=0.;
534 double nu_elem=0.;
535
536 int elem=0;
537 int facei_loc=0,facei=0;
538 int facej_loc=0,facej=0;
539 int dim=0;
540
541 assert(Aij.nb_dim()==3);
542 assert(Aij.dimension(0)==nb_elem_tot);
543 assert(Aij.dimension(1)==nb_faces_elem);
544 assert(Aij.dimension(2)==nb_faces_elem);
545
546 for (elem=0; elem<nb_elem_tot; elem++)
547 {
548 nu_elem=nu(elem);
549 volume=1./volumes(elem);
550 volume*=nu_elem;
551
552 for (facei_loc=0; facei_loc<nb_faces_elem; facei_loc++)
553 {
554 facei=elem_faces(elem,facei_loc);
555
556 signei=1.;
557 if (face_voisins(facei,0)!=elem) signei=-1.;
558
559 for (facej_loc=facei_loc+1; facej_loc<nb_faces_elem; facej_loc++)
560 {
561 facej=elem_faces(elem,facej_loc);
562
563 signej=1.;
564 if (face_voisins(facej,0)!=elem) signej=-1.;
565
566 psc=0.;
567 for (dim=0; dim<Objet_U::dimension; dim++)
568 psc+=face_normales(facei,dim)*face_normales(facej,dim);
569 psc*=signei*signej;
570 psc*=volume;
571
572 Aij(elem,facei_loc,facej_loc)=psc;
573 Aij(elem,facej_loc,facei_loc)=psc;
574 }
575 }
576 }
577}
578
579double Op_Diff_VEF_Face_Stab::calculer_gradients(int facei, const DoubleTab& rij) const
580{
581 const Domaine_VEF& domaine_VEF = le_dom_vef.valeur();
582
583 const int nb_faces_elem=domaine_VEF.domaine().nb_faces_elem();
584
585 const IntTab& elem_faces=domaine_VEF.elem_faces();
586 const IntTab& face_voisins=domaine_VEF.face_voisins();
587 const IntTab& get_num_fac_loc=domaine_VEF.get_num_fac_loc();
588
589 const DoubleTab& face_normales=domaine_VEF.face_normales();
590 const DoubleVect& volumes=domaine_VEF.volumes();
591
592 double volume=0.;
593 double signek=0.;
594 double psc=0.;
595
596 int elem=0;
597 int elem_loc=0;
598 int facei_loc=0;
599 int facek_loc=0,facek=0;
600 int dim=0;
601
602 double mu_ij=0.;
603 for (elem_loc=0; elem_loc<2; elem_loc++)
604 {
605 facei_loc=get_num_fac_loc(facei,elem_loc);
606 elem=face_voisins(facei,elem_loc);
607
608 if (elem!=-1)
609 {
610 volume+=volumes(elem);
611
612 for (facek_loc=0; facek_loc<nb_faces_elem; facek_loc++)
613 if (facek_loc!=facei_loc)
614 {
615 facek=elem_faces(elem,facek_loc);
616
617 signek=1.;
618 if (face_voisins(facek,0)!=elem) signek=-1.;
619
620 psc=0.;
621 for (dim=0; dim<Objet_U::dimension; dim++)
622 psc+=face_normales(facek,dim)*rij(dim);
623 psc*=signek;
624 if (psc<0.) psc*=-1.;
625
626 mu_ij+=psc;
627 }
628 }
629 }
630
631 mu_ij/=volume;
632 mu_ij*=2.;
633 return mu_ij;
634}
635
636void Op_Diff_VEF_Face_Stab::calculer_min(const DoubleTab& inconnueTab, int& dim, DoubleTab& Minima) const
637{
638 const Domaine_VEF& domaine_VEF = le_dom_vef.valeur();
639
640 const int nb_elem_tot=domaine_VEF.nb_elem_tot();
641 const int nb_faces_elem=domaine_VEF.domaine().nb_faces_elem();
642 const int nb_faces_tot=domaine_VEF.nb_faces_tot();
643
644 const IntTab& elem_faces=domaine_VEF.elem_faces();
645
646 int elem=0;
647 int facei_loc=0,facei=0;
648 int facej_loc=0,facej=0;
649
650 int nb_comp = 1;
651 if(inconnueTab.nb_dim()==2) nb_comp=inconnueTab.dimension(1);
652
653 double inc_i=0.;
654 double inc_j=0.;
655
656 const DoubleVect& inconnue= inconnueTab;
657
658 assert(Minima.nb_dim()==1);
659 assert(Minima.dimension(0)==nb_faces_tot);
660 assert(dim<nb_comp);
661
662 for (facei=0; facei<nb_faces_tot; facei++)
663 Minima(facei)=inconnue[facei*nb_comp+dim];
664
665 for (elem=0; elem<nb_elem_tot; elem++)
666 for (facei_loc=0; facei_loc<nb_faces_elem; facei_loc++)
667 {
668 facei=elem_faces(elem,facei_loc);
669
670 inc_i=inconnue[facei*nb_comp+dim];
671 double& mini=Minima(facei);
672
673 for (facej_loc=facei_loc+1; facej_loc<nb_faces_elem; facej_loc++)
674 {
675 facej=elem_faces(elem,facej_loc);
676
677 inc_j=inconnue[facej*nb_comp+dim];
678 double& minj=Minima(facej);
679
680 if (inc_j<mini) mini=inc_j;
681 if (inc_i<minj) minj=inc_i;
682 }
683 }
684
685 const Domaine_Cl_VEF& domaine_Cl_VEF = la_zcl_vef.valeur();
686 const Conds_lim& les_cl = domaine_Cl_VEF.les_conditions_limites();
687
688 const int nb_bords=les_cl.size();
689
690 int num1=0, num2=0;
691 int n_bord=0;
692 int ind_face=0;
693
694 for (n_bord=0; n_bord<nb_bords; n_bord++)
695 {
696 const Cond_lim& la_cl = domaine_Cl_VEF.les_conditions_limites(n_bord);
697 const Front_VF& le_bord = ref_cast(Front_VF,la_cl->frontiere_dis());
698
699 num1 = 0;
700 num2 = le_bord.nb_faces_tot();
701
702 if (sub_type(Periodique,la_cl.valeur()))
703 {
704 const Periodique& la_cl_perio = ref_cast(Periodique,la_cl.valeur());
705
706 for (ind_face=num1; ind_face<num2; ind_face++)
707 {
708 facei=le_bord.num_face(ind_face);
709 facej=le_bord.num_face(la_cl_perio.face_associee(ind_face));
710
711 if (facei<facej)
712 {
713 double& mini=Minima(facei);
714 double& minj=Minima(facej);
715
716 if (mini<minj) minj=mini;
717 if (minj<mini) mini=minj;
718 }
719 }
720 }
721 }
722}
723
724void Op_Diff_VEF_Face_Stab::calculer_max(const DoubleTab& inconnueTab, int& dim, DoubleTab& Maxima) const
725{
726 const Domaine_VEF& domaine_VEF = le_dom_vef.valeur();
727
728 const int nb_elem_tot=domaine_VEF.nb_elem_tot();
729 const int nb_faces_elem=domaine_VEF.domaine().nb_faces_elem();
730 const int nb_faces_tot=domaine_VEF.nb_faces_tot();
731
732 const IntTab& elem_faces=domaine_VEF.elem_faces();
733
734 int elem=0;
735 int facei_loc=0,facei=0;
736 int facej_loc=0,facej=0;
737
738 int nb_comp = 1;
739 if(inconnueTab.nb_dim()==2) nb_comp=inconnueTab.dimension(1);
740
741 double inc_i=0.;
742 double inc_j=0.;
743
744 const DoubleVect& inconnue = inconnueTab;
745
746 assert(Maxima.nb_dim()==1);
747 assert(Maxima.dimension(0)==nb_faces_tot);
748 assert(dim<nb_comp);
749
750 for (facei=0; facei<nb_faces_tot; facei++)
751 Maxima(facei)=inconnue[facei*nb_comp+dim];
752
753 for (elem=0; elem<nb_elem_tot; elem++)
754 for (facei_loc=0; facei_loc<nb_faces_elem; facei_loc++)
755 {
756 facei=elem_faces(elem,facei_loc);
757
758 inc_i=inconnue[facei*nb_comp+dim];
759 double& maxi=Maxima(facei);
760
761 for (facej_loc=facei_loc+1; facej_loc<nb_faces_elem; facej_loc++)
762 {
763 facej=elem_faces(elem,facej_loc);
764
765 inc_j=inconnue[facej*nb_comp+dim];
766 double& maxj=Maxima(facej);
767
768 if (inc_j>maxi) maxi=inc_j;
769 if (inc_i>maxj) maxj=inc_i;
770 }
771 }
772
773 const Domaine_Cl_VEF& domaine_Cl_VEF = la_zcl_vef.valeur();
774 const Conds_lim& les_cl = domaine_Cl_VEF.les_conditions_limites();
775
776 const int nb_bords=les_cl.size();
777
778 int num1=0, num2=0;
779 int n_bord=0;
780 int ind_face=0;
781
782 for (n_bord=0; n_bord<nb_bords; n_bord++)
783 {
784 const Cond_lim& la_cl = domaine_Cl_VEF.les_conditions_limites(n_bord);
785 const Front_VF& le_bord = ref_cast(Front_VF,la_cl->frontiere_dis());
786
787 num1 = 0;
788 num2 = le_bord.nb_faces_tot();
789
790 if (sub_type(Periodique,la_cl.valeur()))
791 {
792 const Periodique& la_cl_perio = ref_cast(Periodique,la_cl.valeur());
793
794 for (ind_face=num1; ind_face<num2; ind_face++)
795 {
796 facei=le_bord.num_face(ind_face);
797 facej=le_bord.num_face(la_cl_perio.face_associee(ind_face));
798
799 if (facei<facej)
800 {
801 double& maxi=Maxima(facei);
802 double& maxj=Maxima(facej);
803
804 if (maxi>maxj) maxj=maxi;
805 if (maxj>maxi) maxi=maxj;
806 }
807 }
808 }
809
810 }
811}
812
814{
816
817 {
818 const Domaine_VEF& domaine_VEF=le_dom_vef.valeur();
819 const Domaine_Cl_VEF& domaine_Cl_VEF=la_zcl_vef.valeur();
820
821 const Conds_lim& les_cl=domaine_Cl_VEF.les_conditions_limites();
822
823 const int nb_bord=les_cl.size();
824 const int nb_faces_tot=domaine_VEF.nb_faces_tot();
825
826 int ind_face=-1;
827
828 is_dirichlet_faces_.resize(nb_faces_tot);
829 is_dirichlet_faces_=0;
830
831 for (int n_bord=0; n_bord<nb_bord; n_bord++)
832 {
833 const Cond_lim& la_cl = domaine_Cl_VEF.les_conditions_limites(n_bord);
834 const Front_VF& le_bord = ref_cast(Front_VF,la_cl->frontiere_dis());
835 int nb_faces_bord_tot=le_bord.nb_faces_tot();
836 int face=-1;
837
838 if ( (sub_type(Dirichlet,la_cl.valeur()))
839 || (sub_type(Dirichlet_homogene,la_cl.valeur()))
840 )
841 for (ind_face=0; ind_face<nb_faces_bord_tot; ind_face++)
842 {
843 face = le_bord.num_face(ind_face);
844 is_dirichlet_faces_(face)=1;
845 }
846 }
847 }
848}
849
850void Op_Diff_VEF_Face_Stab::ajouter_contribution(const DoubleTab& transporte, Matrice_Morse& matrice) const
851{
852 if (!new_jacobian_)
853 Op_Diff_VEF_Face::ajouter_contribution(transporte,matrice);
854 else
855 {
857 // Fill the nu array because matrix assembly with
858 // ajouter_contribution may occur before the first time step
860 const Domaine_Cl_VEF& domaine_Cl_VEF = la_zcl_vef.valeur();
861 const Domaine_VEF& domaine_VEF = le_dom_vef.valeur();
862 const IntTab& elem_faces = domaine_VEF.elem_faces();
863 const IntTab& face_voisins = domaine_VEF.face_voisins();
864
865 int n1 = domaine_VEF.nb_faces();
866 int nb_comp = 1;
867 int nb_dim = transporte.nb_dim();
868
869 DoubleTab nu;
870 int marq=phi_psi_diffuse(equation());
871 const DoubleVect& porosite_elem = equation().milieu().porosite_elem();
872
873 // either div(phi nu grad inco)
874 // or div(nu grad phi inco)
875 // depending on whether phi_psi or psi is diffused
876 modif_par_porosite_si_flag(nu_,nu,!marq,porosite_elem);
877 DoubleVect porosite_eventuelle(equation().milieu().porosite_face());
878 if (!marq)
879 porosite_eventuelle=1;
880
881
882 if(nb_dim==2)
883 nb_comp=transporte.dimension(1);
884
885 int i,j,num_face;
886 int elem1,elem2;
887 int nb_faces_elem = domaine_VEF.domaine().nb_faces_elem();
888 double val;
889
890 int nb_bords=domaine_VEF.nb_front_Cl();
891 for (int n_bord=0; n_bord<nb_bords; n_bord++)
892 {
893 const Cond_lim& la_cl = domaine_Cl_VEF.les_conditions_limites(n_bord);
894 const Front_VF& le_bord = ref_cast(Front_VF,la_cl->frontiere_dis());
895 int num1 = le_bord.num_premiere_face();
896 int num2 = num1 + le_bord.nb_faces();
897
898 if (sub_type(Periodique,la_cl.valeur()))
899 {
900 const Periodique& la_cl_perio = ref_cast(Periodique,la_cl.valeur());
901 int fac_asso;
902 // iterate over only half the periodic faces
903 // the result will be copied to the associated face at the end
904 int num2b=num1+le_bord.nb_faces()/2;
905 for (num_face=num1; num_face<num2b; num_face++)
906 {
907 elem1 = face_voisins(num_face,0);
908 elem2 = face_voisins(num_face,1);
909 fac_asso = la_cl_perio.face_associee(num_face-num1)+num1;
910 for (i=0; i<nb_faces_elem; i++)
911 {
912 if ( (j=elem_faces(elem1,i)) > num_face )
913 {
914 val = viscA(num_face,j,elem1,nu(elem1));
915 if (val<0.) val=0.;
916 for (int nc=0; nc<nb_comp; nc++)
917 {
918 int n0=num_face*nb_comp+nc;
919 int j0=j*nb_comp+nc;
920
921 matrice(n0,n0)+=val*porosite_eventuelle(num_face);
922 matrice(n0,j0)-=val*porosite_eventuelle(j);
923 matrice(j0,n0)-=val*porosite_eventuelle(num_face);
924 matrice(j0,j0)+=val*porosite_eventuelle(j);
925
926 }
927 }
928 if (elem2!=-1)
929 if ( (j=elem_faces(elem2,i)) > num_face )
930 {
931 val= viscA(num_face,j,elem2,nu(elem2));
932 if (val<0.) val=0.;
933 for (int nc=0; nc<nb_comp; nc++)
934 {
935 int n0=num_face*nb_comp+nc;
936 int j0=j*nb_comp+nc;
937 int n0perio=fac_asso*nb_comp+nc;
938 matrice(n0,n0)+=val*porosite_eventuelle(num_face);
939 matrice(n0,j0)-=val*porosite_eventuelle(j);
940 matrice(j0,n0perio)-=val*porosite_eventuelle(num_face);
941 matrice(j0,j0)+=val*porosite_eventuelle(j);
942
943 }
944 }
945 }
946 }
947
948 }
949 else
950 {
951 for (num_face=num1; num_face<num2; num_face++)
952 {
953 elem1 = face_voisins(num_face,0);
954 for (i=0; i<nb_faces_elem; i++)
955 {
956 if ( (j= elem_faces(elem1,i)) > num_face )
957 {
958 val = viscA(num_face,j,elem1,nu(elem1));
959 if (val<0.) val=0.;
960 for (int nc=0; nc<nb_comp; nc++)
961 {
962 int n0=num_face*nb_comp+nc;
963 int j0=j*nb_comp+nc;
964
965 matrice(n0,n0)+=val*porosite_eventuelle(num_face);
966 matrice(n0,j0)-=val*porosite_eventuelle(j);
967 matrice(j0,n0)-=val*porosite_eventuelle(num_face);
968 matrice(j0,j0)+=val*porosite_eventuelle(j);
969
970 }
971 }
972 }
973 }
974 }
975 }
976 int num_premiere_face = domaine_VEF.premiere_face_int();
977 for (num_face=num_premiere_face; num_face<n1; num_face++)
978 {
979 elem1 = face_voisins(num_face,0);
980 elem2 = face_voisins(num_face,1);
981
982 for (i=0; i<nb_faces_elem; i++)
983 {
984 if ( (j=elem_faces(elem1,i)) > num_face )
985 {
986 val = viscA(num_face,j,elem1,nu(elem1));
987 if (val<0.) val=0.;
988 for (int nc=0; nc<nb_comp; nc++)
989 {
990 int n0=num_face*nb_comp+nc;
991 int j0=j*nb_comp+nc;
992
993 matrice(n0,n0)+=val*porosite_eventuelle(num_face);
994 matrice(n0,j0)-=val*porosite_eventuelle(j);
995 matrice(j0,n0)-=val*porosite_eventuelle(num_face);
996 matrice(j0,j0)+=val*porosite_eventuelle(j);
997
998 }
999 }
1000 if (elem2!=-1)
1001 if ( (j=elem_faces(elem2,i)) > num_face )
1002 {
1003 val= viscA(num_face,j,elem2,nu(elem2));
1004 if (val<0.) val=0.;
1005 for (int nc=0; nc<nb_comp; nc++)
1006 {
1007 int n0=num_face*nb_comp+nc;
1008 int j0=j*nb_comp+nc;
1009
1010 matrice(n0,n0)+=val*porosite_eventuelle(num_face);
1011 matrice(n0,j0)-=val*porosite_eventuelle(j);
1012 matrice(j0,n0)-=val*porosite_eventuelle(num_face);
1013 matrice(j0,j0)+=val*porosite_eventuelle(j);
1014
1015 }
1016 }
1017 }
1018 }
1020 }
1021}
1022void Op_Diff_VEF_Face_Stab::ajouter_contribution_multi_scalaire(const DoubleTab& transporte, Matrice_Morse& matrice) const
1023{
1024 if (!new_jacobian_)
1026 else
1027 {
1029 // Fill the nu array because matrix assembly with
1030 // ajouter_contribution may occur before the first time step
1031 remplir_nu(nu_);
1032 const Domaine_Cl_VEF& domaine_Cl_VEF = la_zcl_vef.valeur();
1033 const Domaine_VEF& domaine_VEF = le_dom_vef.valeur();
1034 const IntTab& elem_faces = domaine_VEF.elem_faces();
1035 const IntTab& face_voisins = domaine_VEF.face_voisins();
1036
1037 int n1 = domaine_VEF.nb_faces();
1038 int nb_comp = 1;
1039 int nb_dim = transporte.nb_dim();
1040
1041 DoubleTab nu;
1042 int marq=phi_psi_diffuse(equation());
1043 const DoubleVect& porosite_elem = equation().milieu().porosite_elem();
1044
1045 // either div(phi nu grad inco)
1046 // or div(nu grad phi inco)
1047 // depending on whether phi_psi or psi is diffused
1048 modif_par_porosite_si_flag(nu_,nu,!marq,porosite_elem);
1049 DoubleVect porosite_eventuelle(equation().milieu().porosite_face());
1050 if (!marq)
1051 porosite_eventuelle=1;
1052
1053
1054 if(nb_dim==2)
1055 nb_comp=transporte.dimension(1);
1056
1057 int i,j,num_face;
1058 int elem1,elem2;
1059 int nb_faces_elem = domaine_VEF.domaine().nb_faces_elem();
1060 double val;
1061
1062 int nb_bords=domaine_VEF.nb_front_Cl();
1063 for (int n_bord=0; n_bord<nb_bords; n_bord++)
1064 {
1065 const Cond_lim& la_cl = domaine_Cl_VEF.les_conditions_limites(n_bord);
1066 const Front_VF& le_bord = ref_cast(Front_VF,la_cl->frontiere_dis());
1067 int num1 = le_bord.num_premiere_face();
1068 int num2 = num1 + le_bord.nb_faces();
1069 if (sub_type(Periodique,la_cl.valeur()))
1070 {
1071 const Periodique& la_cl_perio = ref_cast(Periodique,la_cl.valeur());
1072 int fac_asso;
1073 // iterate over only half the periodic faces
1074 // the result will be copied to the associated face at the end
1075 int num2b=num1+le_bord.nb_faces()/2;
1076 for (num_face=num1; num_face<num2b; num_face++)
1077 {
1078 elem1 = face_voisins(num_face,0);
1079 elem2 = face_voisins(num_face,1);
1080 fac_asso = la_cl_perio.face_associee(num_face-num1)+num1;
1081 for (i=0; i<nb_faces_elem; i++)
1082 {
1083 if ( (j=elem_faces(elem1,i)) > num_face )
1084 {
1085 for (int nc=0; nc<nb_comp; nc++)
1086 {
1087 val = viscA(num_face,j,elem1,nu(elem1,nc));
1088 if (val<0.) val=0.;
1089
1090 int n0=num_face*nb_comp+nc;
1091 int j0=j*nb_comp+nc;
1092
1093 matrice(n0,n0)+=val*porosite_eventuelle(num_face);
1094 matrice(n0,j0)-=val*porosite_eventuelle(j);
1095 matrice(j0,n0)-=val*porosite_eventuelle(num_face);
1096 matrice(j0,j0)+=val*porosite_eventuelle(j);
1097
1098 }
1099 }
1100 if (elem2!=-1)
1101 if ( (j=elem_faces(elem2,i)) > num_face )
1102 {
1103 for (int nc=0; nc<nb_comp; nc++)
1104 {
1105 val = viscA(num_face,j,elem2,nu(elem1,nc));
1106 if (val<0.) val=0.;
1107 int n0=num_face*nb_comp+nc;
1108 int j0=j*nb_comp+nc;
1109 int n0perio=fac_asso*nb_comp+nc;
1110 matrice(n0,n0)+=val*porosite_eventuelle(num_face);
1111 matrice(n0,j0)-=val*porosite_eventuelle(j);
1112 matrice(j0,n0perio)-=val*porosite_eventuelle(num_face);
1113 matrice(j0,j0)+=val*porosite_eventuelle(j);
1114
1115 }
1116 }
1117 }
1118 }
1119
1120 }
1121 else
1122 {
1123 for (num_face=num1; num_face<num2; num_face++)
1124 {
1125 elem1 = face_voisins(num_face,0);
1126 for (i=0; i<nb_faces_elem; i++)
1127 {
1128 if ( (j= elem_faces(elem1,i)) > num_face )
1129 {
1130 for (int nc=0; nc<nb_comp; nc++)
1131 {
1132 val = viscA(num_face,j,elem1,nu(elem1,nc));
1133 if (val<0.) val=0.;
1134 int n0=num_face*nb_comp+nc;
1135 int j0=j*nb_comp+nc;
1136
1137 matrice(n0,n0)+=val*porosite_eventuelle(num_face);
1138 matrice(n0,j0)-=val*porosite_eventuelle(j);
1139 matrice(j0,n0)-=val*porosite_eventuelle(num_face);
1140 matrice(j0,j0)+=val*porosite_eventuelle(j);
1141
1142 }
1143 }
1144 }
1145 }
1146 }
1147 }
1148 for (num_face=domaine_VEF.premiere_face_int(); num_face<n1; num_face++)
1149 {
1150 elem1 = face_voisins(num_face,0);
1151 elem2 = face_voisins(num_face,1);
1152
1153 for (i=0; i<nb_faces_elem; i++)
1154 {
1155 if ( (j=elem_faces(elem1,i)) > num_face )
1156 {
1157 for (int nc=0; nc<nb_comp; nc++)
1158 {
1159 val = viscA(num_face,j,elem1,nu(elem1,nc));
1160 if (val<0.) val=0.;
1161 int n0=num_face*nb_comp+nc;
1162 int j0=j*nb_comp+nc;
1163
1164 matrice(n0,n0)+=val*porosite_eventuelle(num_face);
1165 matrice(n0,j0)-=val*porosite_eventuelle(j);
1166 matrice(j0,n0)-=val*porosite_eventuelle(num_face);
1167 matrice(j0,j0)+=val*porosite_eventuelle(j);
1168
1169 }
1170 }
1171 if (elem2!=-1)
1172 if ( (j=elem_faces(elem2,i)) > num_face )
1173 {
1174 for (int nc=0; nc<nb_comp; nc++)
1175 {
1176 val= viscA(num_face,j,elem2,nu(elem2,nc));
1177 if (val<0.) val=0.;
1178 int n0=num_face*nb_comp+nc;
1179 int j0=j*nb_comp+nc;
1180
1181 matrice(n0,n0)+=val*porosite_eventuelle(num_face);
1182 matrice(n0,j0)-=val*porosite_eventuelle(j);
1183 matrice(j0,n0)-=val*porosite_eventuelle(num_face);
1184 matrice(j0,j0)+=val*porosite_eventuelle(j);
1185
1186 }
1187 }
1188 }
1189 }
1191 }
1192}
class Cond_lim Generic class used to represent any class
Definition Cond_lim.h:31
class Conds_lim This class represents a vector of boundary conditions.
Definition Conds_lim.h:32
Classe Dirichlet_homogene This class is the base class of the hierarchy of homogeneous Dirichlet-type...
Dirichlet This class is the base class of the hierarchy of Dirichlet-type boundary conditions.
Definition Dirichlet.h:31
int nb_faces_elem(int=0) const
Returns the number of faces of type i of the geometric elements that make up the domain.
Definition Domaine.h:484
const Cond_lim & les_conditions_limites(int) const
Returns the i-th boundary condition.
class Domaine_VEF
Definition Domaine_VEF.h:53
int nb_faces() const
Returns the total number of faces.
Definition Domaine_VF.h:471
int nb_faces_tot() const
Returns the total number of faces.
Definition Domaine_VF.h:481
virtual double face_normales(int face, int comp) const
Definition Domaine_VF.h:47
double xv(int num_face, int k) const
Definition Domaine_VF.h:76
double volumes(int i) const
Definition Domaine_VF.h:113
const IntTab & get_num_fac_loc() const
Definition Domaine_VF.h:140
virtual double surface(int i) const
Definition Domaine_VF.h:53
int elem_faces(int i, int j) const
Returns the index of the i-th face of element num_elem; the face numbering convention is.
Definition Domaine_VF.h:542
int premiere_face_int() const
A face is internal if and only if it separates two elements.
Definition Domaine_VF.h:463
int face_voisins(int num_face, int i) const
Returns the neighbouring element of num_face in direction i.
Definition Domaine_VF.h:418
int nb_faces_bord() const
Returns the number of faces on which boundary conditions are applied:
Definition Domaine_VF.h:512
int nb_elem_tot() const
int nb_front_Cl() const
const Domaine & domaine() const
Classe Echange_externe_impose: This class represents the special case of the class.
virtual double h_imp(int num) const
Returns the value of the imposed heat exchange coefficient on the i-th component.
virtual double T_ext(int num) const
Returns the value of the imposed temperature on the i-th component of the boundary field.
Class defining operators and methods for all reading operation in an input flow (file,...
Definition Entree.h:42
virtual const Milieu_base & milieu() const =0
class Front_VF
Definition Front_VF.h:36
int nb_faces() const
Definition Front_VF.h:53
int num_premiere_face() const
Definition Front_VF.h:63
int nb_faces_tot() const
Definition Front_VF.h:58
int num_face(const int) const
Definition Front_VF.h:68
Matrice_Morse class - Represents a (sparse) matrix M, not necessarily square,.
DoubleVect & porosite_elem()
Definition Milieu_base.h:58
DoubleVect & porosite_face()
Definition Milieu_base.h:62
const Equation_base & equation() const
Returns the reference to the equation pointed to by MorEqn::mon_equation.
Definition MorEqn.h:62
Classe Neumann_homogene This class is the base class of the hierarchy of homogeneous Neumann-type bou...
Classe Neumann_paroi This boundary condition corresponds to an imposed flux for the.
Neumann_sortie_libre This class represents an open boundary without imposed velocity.
virtual double flux_impose(int i) const
Returns the value of the imposed flux on the i-th component of the field representing the flux at the...
Definition Neumann.cpp:35
static int dimension
Definition Objet_U.h:94
const Nom & que_suis_je() const
Returns the string identifying the class.
Definition Objet_U.cpp:104
virtual Entree & readOn(Entree &)
Reads an Objet_U from an input stream. Virtual method to override.
Definition Objet_U.cpp:289
virtual Sortie & printOn(Sortie &) const
Writes the object to an output stream. Virtual method to override.
Definition Objet_U.cpp:278
class Op_Diff_VEF_Face_Stab
void ajouter_contribution(const DoubleTab &, Matrice_Morse &) const
void calculer_min(const DoubleTab &, int &, DoubleTab &) const
void ajouter_antidiffusion(const DoubleTab &, const DoubleTab &, DoubleTab &) const
void modifie_pour_Cl(const DoubleTab &, DoubleTab &) const
DoubleTab & ajouter(const DoubleTab &, DoubleTab &) const override
void ajouter_diffusion(const DoubleTab &, const DoubleTab &, DoubleTab &) const
void calculer_max(const DoubleTab &, int &, DoubleTab &) const
void calculer_coefficients(const DoubleTab &, DoubleTab &) const
double calculer_gradients(int, const DoubleTab &) const
void ajouter_contribution_multi_scalaire(const DoubleTab &, Matrice_Morse &) const
void ajouter_operateur_centre(const DoubleTab &, const DoubleTab &, DoubleTab &) const
void completer() override
Associates the operator with the domaine_dis, the domaine_Cl_dis, and the unknown of its equation.
class Op_Diff_VEF_Face
void completer() override
Associates the operator with the domaine_dis, the domaine_Cl_dis, and the unknown of its equation.
void ajouter_contribution_multi_scalaire(const DoubleTab &, Matrice_Morse &) const
void ajouter_contribution(const DoubleTab &, Matrice_Morse &) const
int phi_psi_diffuse(const Equation_base &eq) const
Determine whether to compute div(phi nu grad Psi) or div(nu grad Phi psi).
virtual void remplir_nu(DoubleTab &) const
double viscA(int face_i, int face_j, int num_elem, const _TYPE_ &diffu) const
void modifier_matrice_pour_periodique_apres_contribuer(Matrice_Morse &matrice, const Equation_base &) const
Sums the 2 rows of the associated periodic faces, allowing computations in the code without needing t...
void modifier_matrice_pour_periodique_avant_contribuer(Matrice_Morse &matrice, const Equation_base &) const
Divides the coefficients on the periodic face rows by 2 in preparation for applying modifier_matrice_...
void modifier_flux(const Operateur_base &) const
DoubleTab flux_bords_
class Periodique This class represents a periodic boundary condition.
Definition Periodique.h:31
int face_associee(int i) const
Definition Periodique.h:35
static double mp_min(double)
Definition Process.cpp:391
static double mp_max(double)
Definition Process.cpp:379
static void exit(int exit_code=-1)
Exit routine for TRUST within a Kokkos region.
Definition Process.cpp:466
static int je_suis_maitre()
Returns 1 if on the master processor of the current group (i.e. me() == 0), 0 otherwise.
Definition Process.cpp:82
SFichier is to the C++ ofstream class what Sortie is to the C++ ostream class.
Definition SFichier.h:29
Base class for output streams.
Definition Sortie.h:52
virtual void declare_support_masse_volumique(int ok)
The constructor of a derived class that uses the density field must call this function with the value...
Symetrie On symmetry faces, the following properties hold:
Definition Symetrie.h:37
int nb_dim() const
Definition TRUSTTab.h:199
void resize(_SIZE_ n, RESIZE_OPTIONS opt=RESIZE_OPTIONS::COPY_INIT)
Definition TRUSTTab.tpp:469
_SIZE_ dimension(int d) const
Definition TRUSTTab.tpp:133