TRUST 1.9.8
HPC thermohydraulic platform
Loading...
Searching...
No Matches
Op_Diff_VEF_Face.cpp
1/****************************************************************************
2* Copyright (c) 2026, 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.h>
17#include <Champ_P1NC.h>
18#include <Champ_Uniforme.h>
19#include <Periodique.h>
20#include <Symetrie.h>
21#include <Neumann_homogene.h>
22#include <Neumann_paroi.h>
23#include <Echange_externe_impose.h>
24#include <Echange_externe_radiatif.h>
25#include <Neumann_sortie_libre.h>
26#include <Milieu_base.h>
27#include <TRUSTTrav.h>
28#include <Probleme_base.h>
29#include <Navier_Stokes_std.h>
30#include <Porosites_champ.h>
31#include <Device.h>
32#include <Echange_couplage_thermique.h>
33#include <Champ_front_calc_interne.h>
34#include <Robin_VEF.h>
35
36Implemente_instanciable_sans_constructeur(Op_Diff_VEF_Face,"Op_Diff_VEF_P1NC",Op_Diff_VEF_base);
37
42
43Sortie& Op_Diff_VEF_Face::printOn(Sortie& s ) const
44{
45 return s << que_suis_je() ;
46}
47
48Entree& Op_Diff_VEF_Face::readOn(Entree& s )
49{
50 return s ;
51}
52
53/*! @brief associe le champ de diffusivite
54 *
55 */
56void Op_Diff_VEF_Face::associer_diffusivite(const Champ_base& diffu)
57{
58 diffusivite_ = diffu;
59}
60
65
66const Champ_base& Op_Diff_VEF_Face::diffusivite() const
67{
68 return diffusivite_.valeur();
69}
70
71int ma_func_qui_renvoie_int()
72{
73 return 1;
74}
75
76void Op_Diff_VEF_Face::ajouter_cas_scalaire(const DoubleTab& tab_inconnue,
77 DoubleTab& tab_resu, DoubleTab& tab_flux_bords,
78 DoubleTab& tab_nu,
79 const Domaine_Cl_VEF& domaine_Cl_VEF,
80 const Domaine_VEF& domaine_VEF ) const
81{
82 int nb_faces = domaine_VEF.nb_faces();
83 int nb_faces_elem = domaine_VEF.domaine().nb_faces_elem();
84 int nb_bords=domaine_VEF.nb_front_Cl();
85 const int premiere_face_int=domaine_VEF.premiere_face_int();
86
87 {
88 CIntTabView elem_faces = domaine_VEF.elem_faces().view_ro();
89 CIntTabView face_voisins = domaine_VEF.face_voisins().view_ro();
90 CDoubleTabView face_normale = domaine_VEF.face_normales().view_ro();
91 CDoubleArrView inverse_volumes = domaine_VEF.inverse_volumes().view_ro();
92 CDoubleTabView nu = tab_nu.view_ro();
93 CDoubleTabView inconnue = tab_inconnue.view_ro();
94 DoubleArrView flux_bords = static_cast<ArrOfDouble&>(tab_flux_bords).view_rw();
95 DoubleArrView resu = static_cast<ArrOfDouble&>(tab_resu).view_rw();
96 // On traite les faces bord
97 for (int n_bord = 0; n_bord < nb_bords; n_bord++)
98 {
99 const Cond_lim& la_cl = domaine_Cl_VEF.les_conditions_limites(n_bord);
100 const Front_VF& le_bord = ref_cast(Front_VF, la_cl->frontiere_dis());
101 int num1 = 0;
102 int num2 = le_bord.nb_faces_tot();
103 int nb_faces_bord_reel = le_bord.nb_faces();
104 CIntArrView le_bord_num_face = le_bord.num_face().view_ro();
105 if (sub_type(Periodique, la_cl.valeur()))
106 {
107 const Periodique& la_cl_perio = ref_cast(Periodique, la_cl.valeur());
108 CIntArrView face_associee = la_cl_perio.face_associee().view_ro();
109 Kokkos::parallel_for(start_gpu_timer(__KERNEL_NAME__),
110 Kokkos::RangePolicy<>(num1, nb_faces_bord_reel), KOKKOS_LAMBDA(
111 const int ind_face)
112 {
113 int num_face = le_bord_num_face(ind_face);
114 int fac_asso = face_associee(ind_face);
115 fac_asso = le_bord_num_face(fac_asso);
116 for (int kk = 0; kk < 2; kk++)
117 {
118 int elem = face_voisins(num_face, kk);
119 for (int i = 0; i < nb_faces_elem; i++)
120 {
121 int j = elem_faces(elem, i);
122 if (j > num_face && j != fac_asso)
123 {
124 double valA = viscA(num_face, j, elem, nu(elem, 0), face_voisins, face_normale,
125 inverse_volumes);
126 double flux = valA * (inconnue(j, 0) - inconnue(num_face, 0));
127 Kokkos::atomic_add(&resu(num_face), +flux);
128 if (j < nb_faces) // face reelle
129 Kokkos::atomic_add(&resu(j), -0.5 * flux);
130 }
131 }
132 }
133 });
134 end_gpu_timer(__KERNEL_NAME__);
135 }
136 else // Il n'y a qu'une seule composante, donc on traite
137 // une equation scalaire (pas la vitesse) on a pas a utiliser
138 // le tau tangentiel (les lois de paroi thermiques ne calculent pas
139 // d'echange turbulent a la paroi pour l'instant
140 {
141 Kokkos::parallel_for(start_gpu_timer(__KERNEL_NAME__),
142 Kokkos::RangePolicy<>(num1, num2), KOKKOS_LAMBDA(
143 const int ind_face)
144 {
145 int num_face = le_bord_num_face(ind_face);
146 int elem = face_voisins(num_face, 0);
147 for (int i = 0; i < nb_faces_elem; i++)
148 {
149 int j = elem_faces(elem, i);
150 if (j > num_face || num_face >= nb_faces)
151 {
152 double valA = viscA(num_face, j, elem, nu(elem, 0), face_voisins, face_normale,
153 inverse_volumes);
154 double flux = valA * (inconnue(j, 0) - inconnue(num_face, 0));
155 if (num_face < nb_faces) // face reelle
156 {
157 Kokkos::atomic_add(&resu(num_face), +flux);
158 Kokkos::atomic_add(&flux_bords(num_face), -flux);
159 }
160 if (j < nb_faces) // face reelle
161 {
162 Kokkos::atomic_add(&resu(j), -flux);
163 if (j < premiere_face_int)
164 Kokkos::atomic_add(&flux_bords(j), +flux);
165 }
166 }
167 }
168 });
169 end_gpu_timer(__KERNEL_NAME__);
170 }
171 }
172
173 // Faces internes :
174 Kokkos::parallel_for(start_gpu_timer(__KERNEL_NAME__),
175 Kokkos::MDRangePolicy<Kokkos::Rank<2>>({premiere_face_int, 0}, {nb_faces, 2}),
176 KOKKOS_LAMBDA(const int num_face, const int k)
177 {
178 int elem = face_voisins(num_face, k);
179 for (int i = 0; i < nb_faces_elem; i++)
180 {
181 int j = elem_faces(elem, i);
182 if (j > num_face)
183 {
184 int contrib = 1;
185
186 if (j >= nb_faces) // C'est une face virtuelle
187 {
188 int el1 = face_voisins(j, 0);
189 int el2 = face_voisins(j, 1);
190 if ((el1 == -1) || (el2 == -1))
191 contrib = 0;
192 }
193
194 if (contrib)
195 {
196 double valA = viscA(num_face, j, elem, nu(elem, 0), face_voisins,
197 face_normale, inverse_volumes);
198 double flux = valA * (inconnue(j, 0) - inconnue(num_face, 0));
199 Kokkos::atomic_add(&resu(num_face), flux);
200 if (j < nb_faces) // On traite les faces reelles
201 Kokkos::atomic_add(&resu(j), -flux);
202 }
203 }
204 }
205 });
206 end_gpu_timer(__KERNEL_NAME__);
207 }
208
209 // Neumann :
210 for (int n_bord=0; n_bord<nb_bords; n_bord++)
211 {
212 const Cond_lim& la_cl = domaine_Cl_VEF.les_conditions_limites(n_bord);
213 const Front_VF& le_bord = ref_cast(Front_VF,la_cl->frontiere_dis());
214 int ndeb = le_bord.num_premiere_face();
215 int nfin = ndeb + le_bord.nb_faces();
216 if (sub_type(Neumann_paroi,la_cl.valeur()))
217 {
218 const Neumann_paroi& la_cl_paroi = ref_cast(Neumann_paroi, la_cl.valeur());
219 CDoubleArrView surface = domaine_VEF.face_surfaces().view_ro();
220 CDoubleTabView flux_impose = la_cl_paroi.flux_impose().view_ro();
221 DoubleArrView flux_bords = static_cast<ArrOfDouble&>(tab_flux_bords).view_rw();
222 DoubleArrView resu = static_cast<ArrOfDouble&>(tab_resu).view_rw();
223 Kokkos::parallel_for(start_gpu_timer(__KERNEL_NAME__), Kokkos::RangePolicy<>(ndeb, nfin), KOKKOS_LAMBDA(const int face)
224 {
225 double flux = flux_impose(face-ndeb, 0) * surface(face);
226 resu(face) += flux;
227 flux_bords(face) = flux;
228 });
229 end_gpu_timer(__KERNEL_NAME__);
230 }
231 else if (sub_type(Echange_externe_impose,la_cl.valeur()))
232 {
233 const Echange_externe_impose& la_cl_paroi = ref_cast(Echange_externe_impose, la_cl.valeur());
234 const double coeff = COEFF_STEFAN_BOLTZMANN;
235 const bool has_emissivity = la_cl_paroi.has_emissivite();
236 CDoubleArrView surface = domaine_VEF.face_surfaces().view_ro();
237 CDoubleArrView text = static_cast<const ArrOfDouble&>(la_cl_paroi.tab_T_ext()).view_ro();
238 CDoubleArrView himp = static_cast<const ArrOfDouble&>(la_cl_paroi.tab_h_imp()).view_ro();
239 CDoubleArrView eps;
240 if (has_emissivity) eps = static_cast<const ArrOfDouble&>(la_cl_paroi.tab_emissivite()).view_ro();
241 CDoubleArrView inconnue = static_cast<const ArrOfDouble&>(tab_inconnue).view_ro();
242 DoubleArrView resu = static_cast<ArrOfDouble&>(tab_resu).view_rw();
243 DoubleArrView flux_bords = static_cast<ArrOfDouble&>(tab_flux_bords).view_wo();
244 Kokkos::parallel_for(start_gpu_timer(__KERNEL_NAME__), Kokkos::RangePolicy<>(ndeb, nfin), KOKKOS_LAMBDA (const int face)
245 {
246 int ind_face = face - ndeb;
247 double flux = himp(ind_face)*(text(ind_face)-inconnue(face))*surface(face);
248 resu[face] += flux;
249 flux_bords(face) = flux;
250
251 if (has_emissivity)
252 {
253 double T = inconnue(face);
254 double t_ext = text(ind_face);
255 flux = coeff * eps(ind_face) * (t_ext * t_ext * t_ext * t_ext - T * T * T * T) * surface(face);
256 resu[face] += flux;
257 flux_bords(face) += flux;
258 }
259 });
260 end_gpu_timer(__KERNEL_NAME__);
261 }
262 else if (sub_type(Echange_couplage_thermique, la_cl.valeur()))
263 {
264 ToDo_Kokkos("critical");
265 const Echange_couplage_thermique& la_cl_paroi = ref_cast(Echange_couplage_thermique, la_cl.valeur());
266 const DoubleVect& surface = domaine_VEF.face_surfaces();
267 for (int face=ndeb; face<nfin; face++)
268 {
269 double h=la_cl_paroi.h_imp(face-ndeb);
270 double Text=la_cl_paroi.T_ext(face-ndeb);
271 double phiext=la_cl_paroi.flux_exterieur_impose(face-ndeb);
272 double flux=(phiext+h*(Text-tab_inconnue(face)))*surface(face);
273 tab_resu[face] += flux;
274 tab_flux_bords(face) = flux;
275 }
276 }
277 else if (sub_type(Neumann_homogene,la_cl.valeur())
278 || sub_type(Symetrie,la_cl.valeur())
279 || sub_type(Neumann_sortie_libre,la_cl.valeur()))
280 {
281 DoubleArrView flux_bords = static_cast<ArrOfDouble&>(tab_flux_bords).view_wo();
282 Kokkos::parallel_for(start_gpu_timer(__KERNEL_NAME__), Kokkos::RangePolicy<>(ndeb, nfin), KOKKOS_LAMBDA(const int face)
283 {
284 flux_bords(face) = 0.;
285 });
286 end_gpu_timer(__KERNEL_NAME__);
287 }
288 }
289}
290
291void Op_Diff_VEF_Face::ajouter_cas_vectoriel(const DoubleTab& inconnue,
292 DoubleTab& resu, DoubleTab& tab_flux_bords,
293 DoubleTab& nu,
294 const Domaine_Cl_VEF& domaine_Cl_VEF,
295 const Domaine_VEF& domaine_VEF,
296 int nb_comp) const
297{
298 assert(nb_comp==dimension);
299
300 // Construction du tableau grad_ si necessaire
301 if(!grad_.get_md_vector())
302 {
303 grad_.resize(0, Objet_U::dimension, Objet_U::dimension);
304 domaine_VEF.domaine().creer_tableau_elements(grad_);
305 }
306 Champ_P1NC::calcul_gradient(inconnue,grad_,domaine_Cl_VEF);
307
308 /* ToDo OpenMP : factoriser avec Op_Dift_VEF_Face.cpp dans une classe template
309 if (le_modele_turbulence->utiliser_loi_paroi())
310 {
311 Champ_P1NC::calcul_duidxj_paroi(grad_,nu,nu_turb,tau_tan_,domaine_Cl_VEF);
312 grad_.echange_espace_virtuel(); // gradient_elem a jour sur les elements virtuels
313 }
314 DoubleTab Re;
315 Re.resize(0, Objet_U::dimension, Objet_U::dimension);
316 domaine_VEF.domaine().creer_tableau_elements(Re);
317 Re = 0.;
318 if (le_modele_turbulence->calcul_tenseur_Re(nu_turb, grad_, Re))
319 {
320 Cerr << "On utilise une diffusion turbulente non linaire dans NS" << finl;
321 for (int elem=0; elem<nb_elem; elem++)
322 for (int i=0; i<nbr_comp; i++)
323 for (int j=0; j<nbr_comp; j++)
324 Re(elem,i,j) *= nu_turb[elem];
325 }
326 else
327 {
328 for (int elem=0; elem<nb_elem; elem++)
329 for (int i=0; i<nbr_comp; i++)
330 for (int j=0; j<nbr_comp; j++)
331 Re(elem,i,j) = nu_turb[elem]*(grad_(elem,i,j) + grad_(elem,j,i));
332 }
333 Re.echange_espace_virtuel();
334 */
335
336 int nb_faces = domaine_VEF.nb_faces();
337 int nb_faces_bord = domaine_VEF.premiere_face_int();
338 CIntTabView face_voisins_v = domaine_VEF.face_voisins().view_ro();
339 CDoubleTabView face_normales_v = domaine_VEF.face_normales().view_ro();
340 CDoubleTabView nu_v = nu.view_ro();
341 CDoubleTabView3 grad_v = grad_.view_ro<3>();
342 DoubleTabView resu_v = resu.view_rw();
343 DoubleTabView tab_flux_bords_v = tab_flux_bords.view_rw();
344
345 auto kern_ajouter = KOKKOS_LAMBDA(int
346 num_face, int k)
347 {
348 int elem = face_voisins_v(num_face, k);
349 if (elem >= 0)
350 {
351 int ori = 1 - 2 * k;
352 double nu_elem = nu_v(elem, 0);
353 for (int i = 0; i < nb_comp; i++)
354 for (int j = 0; j < nb_comp; j++)
355 {
356 double grad_ij = grad_v(elem, i, j);
357 double grad_ji = grad_v(elem, j, i);
358 double fn = face_normales_v(num_face, j);
359 double flux = ori * fn * (nu_elem * grad_ij /* + Re(elem, i, j) */ );
360 Kokkos::atomic_sub(&resu_v(num_face, i), flux);
361
362 if (num_face < nb_faces_bord)
363 {
364 double flux_bord = ori * fn * (nu_elem * (grad_ij + grad_ji));
365 Kokkos::atomic_sub(&tab_flux_bords_v(num_face, i), flux_bord);
366 }
367 }
368 }
369 };
370 Kokkos::parallel_for(start_gpu_timer(__KERNEL_NAME__), Kokkos::MDRangePolicy<Kokkos::Rank<2>>({0,0}, {nb_faces,2}) , kern_ajouter);
371 end_gpu_timer(__KERNEL_NAME__);
372
373
374 const int nb_bords=domaine_VEF.nb_front_Cl();
375 for (int n_bord=0; n_bord<nb_bords; n_bord++)
376 {
377 // Update flux_bords on symmetry:
378 const Cond_lim& la_cl = domaine_Cl_VEF.les_conditions_limites(n_bord);
379 if (sub_type(Symetrie,la_cl.valeur()))
380 {
381 const Front_VF& le_bord = ref_cast(Front_VF,la_cl->frontiere_dis());
382 int ndeb = le_bord.num_premiere_face();
383 int nfin = ndeb + le_bord.nb_faces();
384 DoubleTabView flux_bords = tab_flux_bords.view_wo();
385 Kokkos::parallel_for(start_gpu_timer(__KERNEL_NAME__), Kokkos::RangePolicy<>(ndeb, nfin), KOKKOS_LAMBDA(const int face)
386 {
387 flux_bords(face, 0) = 0.;
388 });
389 end_gpu_timer(__KERNEL_NAME__);
390 }
391
392 else if (sub_type(Robin_VEF, la_cl.valeur()))
393 {
394#ifdef TRUST_USE_GPU
395 Cerr << "Warning not tested on GPU" << finl;
396#endif
397 ToDo_Kokkos("critical");
398
399 const Front_VF& le_bord = ref_cast(Front_VF,la_cl->frontiere_dis());
400 const Robin_VEF& la_cl_robin = ref_cast(Robin_VEF,la_cl.valeur());
401 int marq = phi_psi_diffuse(equation());
402 const DoubleVect& porosite_face = equation().milieu().porosite_face();
403 double scale_factor_is_one = 1.;
404 double inv_alpha = 1./la_cl_robin.get_alpha_cl() ;
405 double inv_beta = 1./la_cl_robin.get_beta_cl();
406 double inv_alpha_minus_inv_beta = 1./la_cl_robin.get_alpha_cl() - 1./la_cl_robin.get_beta_cl();
407 DoubleTab normal_vector;
408 int ndeb = le_bord.num_premiere_face();
409 int nfin = ndeb +le_bord.nb_faces();
410 for (int face=ndeb; face<nfin; face++)
411 {
412 int id_face_bord = face -ndeb;
413 double face_surface = domaine_VEF.face_surfaces(face);
414 normal_vector = domaine_VEF.normalized_boundaries_outward_vector(face, scale_factor_is_one);
415
416 for (int nc1=0; nc1<nb_comp; nc1++)
417 {
418 double flux_robin_uu = 0. ;
419 double flux_robin_rhs = 0.;
420 double flux_tot ;
421
422 // forward term for velocity
423 for (int nc2 = 0; nc2<nb_comp; nc2++)
424 {
425 const double normal2 = normal_vector(nc1)*normal_vector(nc2);
426 flux_robin_uu += (inv_beta*(nc1==nc2) + inv_alpha_minus_inv_beta*normal2)* (face_surface);
427 }
428 flux_robin_uu *= inconnue(face,nc1);
429
430 // rhs for robin bc
431 double val;
432 if (dimension == 2)
433 {
434 // add normal component rhs
435 val = inv_alpha * normal_vector(nc1) * la_cl_robin.flux_normal_imp(id_face_bord);
436
437 // add tangential component rhs
438 double tgte = (2*nc1-1)*normal_vector(1-nc1);
439 val += inv_beta * la_cl_robin.flux_tangentiel_imp(id_face_bord, 0)*tgte;
440 }
441 else
442 {
443 // add normal component rhs
444 val = inv_alpha * normal_vector(nc1) * la_cl_robin.flux_normal_imp(id_face_bord);
445
446 // add tangential component rhs
447 val += inv_beta * la_cl_robin.flux_tangentiel_imp(id_face_bord, nc1) ;
448 }
449 flux_robin_rhs = val*face_surface;
450 flux_tot = (flux_robin_rhs - flux_robin_uu)* (marq ? porosite_face(face) : 1);
451 resu(face,nc1) += flux_tot;
452 tab_flux_bords(face,nc1) += flux_tot;
453
454
455 }
456 }
457 }
458 }
459}
460
461void Op_Diff_VEF_Face::ajouter_cas_multi_scalaire(const DoubleTab& inconnue,
462 DoubleTab& resu, DoubleTab& tab_flux_bords,
463 DoubleTab& nu,
464 const Domaine_Cl_VEF& domaine_Cl_VEF,
465 const Domaine_VEF& domaine_VEF,
466 int nb_comp) const
467{
468 ToDo_Kokkos("critical");
469 const IntTab& elemfaces = domaine_VEF.elem_faces();
470 const IntTab& face_voisins = domaine_VEF.face_voisins();
471 int i0,j,num_face;
472 int nb_faces = domaine_VEF.nb_faces();
473 int nb_faces_elem = domaine_VEF.domaine().nb_faces_elem();
474 int n_bord;
475 double flux0;
476 //DoubleVect n(Objet_U::dimension);
477 //DoubleTrav Tgrad(Objet_U::dimension,Objet_U::dimension);
478
479 assert(nb_comp>1);
480 int nb_bords=domaine_VEF.nb_front_Cl();
481 int ind_face;
482
483 for (n_bord=0; n_bord<nb_bords; n_bord++)
484 {
485 const Cond_lim& la_cl = domaine_Cl_VEF.les_conditions_limites(n_bord);
486 const Front_VF& le_bord = ref_cast(Front_VF,la_cl->frontiere_dis());
487 // const IntTab& elemfaces = domaine_VEF.elem_faces();
488 int num1 = 0;
489 int num2 = le_bord.nb_faces_tot();
490 int nb_faces_bord_reel = le_bord.nb_faces();
491
492 if (sub_type(Periodique,la_cl.valeur()))
493 {
494 const Periodique& la_cl_perio = ref_cast(Periodique,la_cl.valeur());
495 int fac_asso;
496 for (ind_face=num1; ind_face<nb_faces_bord_reel; ind_face++)
497 {
498 fac_asso = la_cl_perio.face_associee(ind_face);
499 fac_asso = le_bord.num_face(fac_asso);
500 num_face = le_bord.num_face(ind_face);
501 for (int kk=0; kk<2; kk++)
502 {
503 int elem = face_voisins(num_face, kk);
504 for (i0=0; i0<nb_faces_elem; i0++)
505 {
506 if ( ( (j= elemfaces(elem,i0)) > num_face ) && (j != fac_asso ) )
507 {
508 for (int nc=0; nc<nb_comp; nc++)
509 {
510 double valA = viscA(num_face,j,elem,nu(elem,nc));
511 resu(num_face,nc)+=valA*inconnue(j,nc);
512 resu(num_face,nc)-=valA*inconnue(num_face,nc);
513 if(j<nb_faces) // face reelle
514 {
515 ////ATENTION DIFF NUM_face avec ma version
516 resu(j,nc)+=0.5*valA*inconnue(num_face,nc);
517 resu(j,nc)-=0.5*valA*inconnue(j,nc);
518 }
519 }
520 }
521 }
522 }
523 }
524 }// fin if periodique
525 else
526 {
527 for (ind_face=num1; ind_face<num2; ind_face++)
528 {
529 num_face = le_bord.num_face(ind_face);
530 int elem=face_voisins(num_face,0);
531
532 // Boucle sur les faces :
533 for (int i=0; i<nb_faces_elem; i++)
534 if (( (j= elemfaces(elem,i)) > num_face ) || (ind_face>=nb_faces_bord_reel))
535 {
536 for (int nc=0; nc<nb_comp; nc++)
537 {
538 double valA = viscA(num_face,j,elem,nu(elem,nc));
539 if (ind_face<nb_faces_bord_reel)
540 {
541 double flux=valA*(inconnue(j,nc)-inconnue(num_face,nc));
542 resu(num_face,nc)+=flux;
543 tab_flux_bords(num_face,nc)-=flux;
544 }
545
546 if(j<nb_faces) // face reelle
547 {
548 resu(j,nc)+=valA*inconnue(num_face,nc);
549 resu(j,nc)-=valA*inconnue(j,nc);
550 }
551 }
552 }
553 }
554 }
555 }//Fin for n_bord
556
557 // On traite les faces internes
558
559 for (num_face=domaine_VEF.premiere_face_int(); num_face<nb_faces; num_face++)
560 {
561 for (int k=0; k<2; k++)
562 {
563 int elem = face_voisins(num_face,k);
564 for (i0=0; i0<nb_faces_elem; i0++)
565 {
566 if ( (j= elemfaces(elem,i0)) > num_face )
567 {
568 int el1,el2;
569 int contrib=1;
570 if(j>=nb_faces) // C'est une face virtuelle
571 {
572 el1 = face_voisins(j,0);
573 el2 = face_voisins(j,1);
574 if((el1==-1)||(el2==-1))
575 contrib=0;
576 }
577 if(contrib)
578 {
579 for (int nc=0; nc<nb_comp; nc++)
580 {
581 double valA = viscA(num_face,j,elem,nu(elem,nc));
582 resu(num_face,nc)+=valA*inconnue(j,nc);
583 resu(num_face,nc)-=valA*inconnue(num_face,nc);
584 if(j<nb_faces) // On traite les faces reelles
585 {
586 resu(j,nc)+=valA*inconnue(num_face,nc);
587 resu(j,nc)-=valA*inconnue(j,nc);
588 }
589 else
590 {
591 // La face j est virtuelle
592 }
593 }
594 }
595 }
596 }
597 }
598 }// Fin faces internes
599
600
601 //On se base sur ce qui est fait pour le cas scalaire
602 for (n_bord=0; n_bord<nb_bords; n_bord++)
603 {
604 const Cond_lim& la_cl = domaine_Cl_VEF.les_conditions_limites(n_bord);
605
606 if (sub_type(Neumann_paroi,la_cl.valeur()))
607 {
608 const Neumann_paroi& la_cl_paroi = ref_cast(Neumann_paroi, la_cl.valeur());
609 const Front_VF& le_bord = ref_cast(Front_VF,la_cl->frontiere_dis());
610 int ndeb = le_bord.num_premiere_face();
611 int nfin = ndeb + le_bord.nb_faces();
612 for (int face=ndeb; face<nfin; face++)
613 {
614 for (int nc=0; nc<nb_comp; nc++)
615 {
616 flux0=la_cl_paroi.flux_impose(face-ndeb,nc)*domaine_VEF.surface(face);
617 resu(face,nc) += flux0;
618 tab_flux_bords(face,nc) = flux0;
619 }
620 }
621 }
622 else if (sub_type(Echange_externe_impose,la_cl.valeur()))
623 {
624 throw;
625 const Echange_externe_impose& la_cl_paroi = ref_cast(Echange_externe_impose, la_cl.valeur());
626 const Front_VF& le_bord = ref_cast(Front_VF,la_cl->frontiere_dis());
627 int ndeb = le_bord.num_premiere_face();
628 int nfin = ndeb + le_bord.nb_faces();
629 for (int face=ndeb; face<nfin; face++)
630 {
631 for (int nc=0; nc<nb_comp; nc++)
632 {
633 flux0=la_cl_paroi.h_imp(face-ndeb,nc)*(la_cl_paroi.T_ext(face-ndeb,nc)-inconnue(face,nc))*domaine_VEF.surface(face);
634 resu(face,nc) += flux0;
635 tab_flux_bords(face,nc) = flux0;
636
637 if (la_cl_paroi.has_emissivite())
638 {
639 const double text = la_cl_paroi.T_ext(face - ndeb, nc), T = inconnue(face, nc);
640 flux0 = COEFF_STEFAN_BOLTZMANN * la_cl_paroi.emissivite(face - ndeb, nc) * (text * text * text * text - T * T * T * T) * domaine_VEF.face_surfaces(face);
641 resu(face, nc) += flux0;
642 tab_flux_bords(face, nc) += flux0;
643 }
644 }
645 }
646 }
647 else if (sub_type(Echange_couplage_thermique,la_cl.valeur()))
648 {
649 const Echange_couplage_thermique& la_cl_paroi = ref_cast(Echange_couplage_thermique, la_cl.valeur());
650 const Front_VF& le_bord = ref_cast(Front_VF,la_cl->frontiere_dis());
651 int ndeb = le_bord.num_premiere_face();
652 int nfin = ndeb + le_bord.nb_faces();
653 for (int face=ndeb; face<nfin; face++)
654 {
655 for (int nc=0; nc<nb_comp; nc++)
656 {
657 double phiext = la_cl_paroi.flux_exterieur_impose(face-ndeb,nc);
658 flux0 = (phiext + la_cl_paroi.h_imp(face-ndeb,nc)*(la_cl_paroi.T_ext(face-ndeb,nc)-inconnue(face,nc)))*domaine_VEF.surface(face);
659 resu(face,nc) += flux0;
660 tab_flux_bords(face,nc) = flux0;
661 }
662 }
663 }
664 else if (sub_type(Neumann_homogene,la_cl.valeur())
665 || sub_type(Symetrie,la_cl.valeur())
666 || sub_type(Neumann_sortie_libre,la_cl.valeur()))
667 {
668 const Front_VF& le_bord = ref_cast(Front_VF,la_cl->frontiere_dis());
669 int ndeb = le_bord.num_premiere_face();
670 int nfin = ndeb + le_bord.nb_faces();
671 for (int face=ndeb; face<nfin; face++)
672 for (int nc=0; nc<nb_comp; nc++)
673 tab_flux_bords(face,nc) = 0.;
674 }
675 }
676}
677
678
679DoubleTab& Op_Diff_VEF_Face::ajouter(const DoubleTab& inconnue_org, DoubleTab& resu) const
680{
681 remplir_nu(nu_);
682 const Domaine_Cl_VEF& domaine_Cl_VEF = la_zcl_vef.valeur();
683 const Domaine_VEF& domaine_VEF = le_dom_vef.valeur();
684
685 int nb_comp = 1;
686 int nb_dim = resu.nb_dim();
687 if(nb_dim==2)
688 nb_comp=resu.dimension(1);
689 DoubleTab nu;
690 DoubleTab tab_inconnue;
691 int marq=phi_psi_diffuse(equation());
692 const DoubleVect& porosite_face = equation().milieu().porosite_face();
693 const DoubleVect& porosite_elem = equation().milieu().porosite_elem();
694 // soit on a div(phi nu grad inco)
695 // soit on a div(nu grad phi inco)
696 // cela depend si on diffuse phi_psi ou psi
697 modif_par_porosite_si_flag(nu_,nu,!marq,porosite_elem);
698 const DoubleTab& inconnue=modif_par_porosite_si_flag(inconnue_org,tab_inconnue,marq,porosite_face);
699
700 const Champ_base& inco = equation().inconnue();
701 const Nature_du_champ nature_champ = inco.nature_du_champ();
702
703 // On dimensionne et initialise le tableau des bilans de flux:
704 if (flux_bords_.size_array()!=domaine_VEF.nb_faces_bord()) flux_bords_.resize(domaine_VEF.nb_faces_bord(),nature_champ==scalaire ? 1 : nb_comp);
705 flux_bords_=0.;
706
707 if(nature_champ==scalaire)
708 ajouter_cas_scalaire(inconnue, resu, flux_bords_, nu, domaine_Cl_VEF, domaine_VEF);
709 else if (nature_champ==vectoriel)
710 ajouter_cas_vectoriel(inconnue, resu, flux_bords_, nu, domaine_Cl_VEF, domaine_VEF,nb_comp);
711 else if (nature_champ==multi_scalaire)
712 ajouter_cas_multi_scalaire(inconnue, resu, flux_bords_, nu, domaine_Cl_VEF, domaine_VEF,nb_comp);
713 modifier_flux(*this);
714
715 return resu;
716}
717
718DoubleTab& Op_Diff_VEF_Face::calculer(const DoubleTab& inconnue, DoubleTab& resu) const
719{
720 resu = 0;
721 return ajouter(inconnue,resu);
722}
723
724void Op_Diff_VEF_Face::ajouter_contribution(const DoubleTab& tab_transporte, Matrice_Morse& tab_matrice) const
725{
726 const Domaine_Cl_VEF& domaine_Cl_VEF = la_zcl_vef.valeur();
727 const Domaine_VEF& domaine_VEF = le_dom_vef.valeur();
728
729 modifier_matrice_pour_periodique_avant_contribuer(tab_matrice,equation());
730
731 // On remplit le tableau nu car l'assemblage d'une
732 // matrice avec ajouter_contribution peut se faire
733 // avant le premier pas de temps
734 remplir_nu(nu_);
735 DoubleTrav tab_nu;
736
737 // soit on a div(phi nu grad inco)
738 // soit on a div(nu grad phi inco)
739 // cela depend si on diffuse phi_psi ou psi
740 int marq = phi_psi_diffuse(equation());
741 modif_par_porosite_si_flag(nu_,tab_nu,!marq,equation().milieu().porosite_elem());
742
743 int nb_dim = tab_transporte.nb_dim();
744 int nb_comp = (nb_dim==2 ? tab_transporte.dimension(1) : 1);
745 int nb_faces_elem = domaine_VEF.domaine().nb_faces_elem();
746 int nb_bords = domaine_VEF.nb_front_Cl();
747
748
749 IntTrav tab_face_associee(domaine_VEF.premiere_face_int());
750 IntTrav tab_fac2b_idx(domaine_VEF.nb_faces());
751
752
753 // Recuperer les indices des faces de bord periodiques et les
754 // faces associees dans des tableaux au prealable pour
755 // permettre un acces structure dans le kernel ensuite
756 for (int n_bord = 0; n_bord < nb_bords; n_bord++)
757 {
758 const Cond_lim& la_cl = domaine_Cl_VEF.les_conditions_limites(n_bord);
759 const Front_VF& le_bord = ref_cast(Front_VF,la_cl->frontiere_dis());
760 int num1 = le_bord.num_premiere_face();
761
762 if (sub_type(Periodique, la_cl.valeur()))
763 {
764 const Periodique& la_cl_perio = ref_cast(Periodique, la_cl.valeur());
765
766 int nb_faces = le_bord.nb_faces();
767 int num2b = num1 + nb_faces / 2;
768
769 // on ne parcourt que la moitie des faces periodiques
770 // on copiera a la fin le resultat dans la face
771 // associee...
772 ToDo_Kokkos("critical");
773 for (int fac = num1; fac < num2b; fac++)
774 {
775 int fac_asso = la_cl_perio.face_associee(fac - num1) + num1;
776
777 tab_face_associee(fac) = fac_asso;
778 tab_fac2b_idx(fac) = fac;
779 tab_fac2b_idx(fac+nb_faces/2) = fac;
780 }
781 }
782 }
783
784 CIntTabView elem_faces = domaine_VEF.elem_faces().view_ro();
785 CIntTabView face_voisins = domaine_VEF.face_voisins().view_ro();
786 CDoubleArrView porosite_face = equation().milieu().porosite_face().view_ro();
787 CDoubleArrView inverse_volumes = domaine_VEF.inverse_volumes().view_ro();
788 CDoubleTabView face_normale = domaine_VEF.face_normales().view_ro();
789 CDoubleTabView nu = tab_nu.view_ro();
790 CIntArrView est_face_bord = domaine_VEF.est_face_bord().view_ro();
791 CIntArrView face_associee = static_cast<ArrOfInt&>(tab_face_associee).view_ro();
792 CIntArrView fac2b_idx = static_cast<ArrOfInt&>(tab_fac2b_idx).view_ro();
793 Matrice_Morse_View matrice;
794 matrice.set(tab_matrice);
795
796 auto ajouter_contrib = KOKKOS_LAMBDA(const int fac)
797 {
798 int type_face = est_face_bord(fac);
799 int fac2b = fac2b_idx(fac);
800
801 if (type_face == 2 && fac == fac2b) // faces perio
802 {
803 int elem1 = face_voisins(fac,0);
804 int elem2 = face_voisins(fac,1);
805
806 int fac_asso = face_associee(fac);
807 for (int i = 0; i < nb_faces_elem; i++)
808 {
809 int j = elem_faces(elem1,i);
810 if (j > fac)
811 {
812 double val = viscA(fac,j,elem1,nu(elem1,0), face_voisins, face_normale, inverse_volumes);
813 double coeff_face1 = val * (marq ? porosite_face(fac) : 1);
814 double coeff_face2 = val * (marq ? porosite_face(j) : 1);
815
816 for (int nc = 0; nc < nb_comp; nc++)
817 {
818 int n0 = fac*nb_comp + nc;
819 int j0 = j*nb_comp + nc;
820
821 matrice.atomic_add(n0, n0, +coeff_face1);
822 matrice.atomic_add(n0, j0, -coeff_face2);
823 matrice.atomic_add(j0, n0, -coeff_face1);
824 matrice.atomic_add(j0, j0, +coeff_face2);
825 }
826 }
827 if (elem2 != -1)
828 {
829 j = elem_faces(elem2,i);
830 if (j > fac)
831 {
832 double val = viscA(fac,j,elem2,nu(elem2,0), face_voisins, face_normale, inverse_volumes);
833 double coeff_face1 = val * (marq ? porosite_face(fac) : 1);
834 double coeff_face2 = val * (marq ? porosite_face(j) : 1);
835
836 for (int nc = 0; nc < nb_comp; nc++)
837 {
838 int n0 = fac*nb_comp + nc;
839 int j0 = j*nb_comp + nc;
840 int n1 = fac_asso*nb_comp+nc;
841
842 matrice.atomic_add(n0, n0, +coeff_face1);
843 matrice.atomic_add(n0, j0, -coeff_face2);
844 matrice.atomic_add(j0, n1, -coeff_face1);
845 matrice.atomic_add(j0, j0, +coeff_face2);
846 }
847 }
848 }
849 }
850 }
851 else if (type_face == 1) // faces bord non perio
852 {
853 int elem1 = face_voisins(fac,0);
854 for (int i = 0; i < nb_faces_elem; i++)
855 {
856 int j = elem_faces(elem1,i);
857 if (j > fac)
858 {
859 double val = viscA(fac,j,elem1,nu(elem1,0), face_voisins, face_normale, inverse_volumes);
860 double coeff_face1 = val * (marq ? porosite_face(fac) : 1);
861 double coeff_face2 = val * (marq ? porosite_face(j) : 1);
862
863 for (int nc = 0; nc < nb_comp; nc++)
864 {
865 int n0 = fac*nb_comp + nc;
866 int j0 = j*nb_comp + nc;
867
868 matrice.atomic_add(n0, n0, +coeff_face1);
869 matrice.atomic_add(n0, j0, -coeff_face2);
870 matrice.atomic_add(j0, n0, -coeff_face1);
871 matrice.atomic_add(j0, j0, +coeff_face2);
872 }
873 }
874 }
875 }
876 else if (type_face == 0) // faces internes
877 {
878 for (int k=0; k<2; k++)
879 {
880 int elem = face_voisins(fac,k);
881 if (elem!=-1)
882 {
883 for (int i = 0; i < nb_faces_elem; i++)
884 {
885 int j = elem_faces(elem, i);
886 if (j > fac)
887 {
888 double val = viscA(fac, j, elem, nu(elem, 0), face_voisins, face_normale, inverse_volumes);
889 double coeff_face1 = val * (marq ? porosite_face(fac) : 1);
890 double coeff_face2 = val * (marq ? porosite_face(j) : 1);
891
892 for (int nc = 0; nc < nb_comp; nc++)
893 {
894 int n0 = fac * nb_comp + nc;
895 int j0 = j * nb_comp + nc;
896
897 matrice.atomic_add(n0, n0, +coeff_face1);
898 matrice.atomic_add(n0, j0, -coeff_face2);
899 matrice.atomic_add(j0, n0, -coeff_face1);
900 matrice.atomic_add(j0, j0, +coeff_face2);
901 }
902 }
903 }
904 }
905 }
906 }
907 };
908 Kokkos::parallel_for(start_gpu_timer(__KERNEL_NAME__),domaine_VEF.nb_faces(), ajouter_contrib);
909 end_gpu_timer(__KERNEL_NAME__);
910
911 int premiere_face_int = domaine_VEF.premiere_face_int();
912 DoubleTrav tab_h_impose(premiere_face_int);
913 DoubleTrav tab_derivee_flux_exterieur_imposee(premiere_face_int);
914
915 // Neumann: remplir les tableaux avec les conditions aux
916 // bords pour le kernel Kokkos
917 for (int n_bord = 0; n_bord < nb_bords; n_bord++)
918 {
919 const Cond_lim& la_cl = domaine_Cl_VEF.les_conditions_limites(n_bord);
920
921 if (sub_type(Echange_externe_impose,la_cl.valeur()))
922 {
923 const Echange_externe_impose& la_cl_paroi = ref_cast(Echange_externe_impose, la_cl.valeur());
924 const Front_VF& le_bord = ref_cast(Front_VF,la_cl->frontiere_dis());
925 int ndeb = le_bord.num_premiere_face();
926 int nfin = ndeb + le_bord.nb_faces();
927 const double coeff = COEFF_STEFAN_BOLTZMANN;
928 const bool has_emissivity = la_cl_paroi.has_emissivite();
929 CDoubleArrView inconnue = static_cast<const ArrOfDouble&>(equation().inconnue().valeurs()).view_ro();
930 CDoubleArrView himp = static_cast<const ArrOfDouble&>(la_cl_paroi.tab_h_imp()).view_ro();
931 CDoubleArrView eps;
932 if (has_emissivity) eps = static_cast<const ArrOfDouble&>(la_cl_paroi.tab_emissivite()).view_ro();
933 DoubleArrView h_impose = static_cast<ArrOfDouble&>(tab_h_impose).view_wo();
934 Kokkos::parallel_for(start_gpu_timer(__KERNEL_NAME__), Kokkos::RangePolicy<>(ndeb, nfin), KOKKOS_LAMBDA(const int face)
935 {
936 int ind_face = face - ndeb;
937 h_impose(face) = himp(ind_face);
938 if (has_emissivity)
939 {
940 const double T = inconnue(face);
941 h_impose(face) = 4 * coeff * eps(ind_face) * T * T * T;
942 }
943 });
944 end_gpu_timer(__KERNEL_NAME__);
945 }
946 else if (sub_type(Echange_couplage_thermique, la_cl.valeur()))
947 {
948 const Echange_couplage_thermique& la_cl_paroi = ref_cast(Echange_couplage_thermique, la_cl.valeur());
949 const Front_VF& le_bord = ref_cast(Front_VF,la_cl->frontiere_dis());
950
951 int ndeb = le_bord.num_premiere_face();
952 int nfin = ndeb + le_bord.nb_faces();
953 ToDo_Kokkos("critical");
954 for (int face = ndeb; face < nfin; face++)
955 {
956 tab_h_impose(face) = la_cl_paroi.h_imp(face-ndeb);
957 tab_derivee_flux_exterieur_imposee(face) = la_cl_paroi.derivee_flux_exterieur_imposee(face-ndeb);
958 }
959 }
960 else if (sub_type(Robin_VEF, la_cl.valeur()))
961 {
962 ToDo_Kokkos("critical");
963 const Front_VF& le_bord = ref_cast(Front_VF,la_cl->frontiere_dis());
964 const Robin_VEF& la_cl_robin = ref_cast(Robin_VEF,la_cl.valeur());
965 double inv_alpha_minus_inv_beta = 1./la_cl_robin.get_alpha_cl() - 1./la_cl_robin.get_beta_cl();
966 double inv_beta = 1./la_cl_robin.get_beta_cl();
967 double scale_factor_is_one = 1. ;
968 DoubleTab normal_vector ;
969 int ndeb = le_bord.num_premiere_face();
970 int nfin = ndeb + le_bord.nb_faces();
971 const DoubleVect& tab_porosite_face = equation().milieu().porosite_face();
972 for (int face = ndeb; face < nfin; face++)
973 {
974 double face_surface = domaine_VEF.face_surfaces(face);
975 normal_vector = domaine_VEF.normalized_boundaries_outward_vector(face, scale_factor_is_one);
976 //int elem = face_voisins(face, 0) ;
977 for (int nc1 = 0; nc1 < nb_comp; nc1++)
978 {
979 const int i = face * nb_comp + nc1;
980
981 // diagonal term
982 double val = (inv_beta + inv_alpha_minus_inv_beta*normal_vector(nc1)*normal_vector(nc1))*face_surface;
983 double robin_contribution= val * (marq ? tab_porosite_face(face) : 1) ;
984 tab_matrice(i,i) += robin_contribution ;
985
986 // extradiagonal term
987 for (int nc2 = 0; nc2<nc1; nc2++)
988 {
989 const int j = face * nb_comp + nc2;
990 const double normal2 = normal_vector(nc1)*normal_vector(nc2);
991 val = inv_alpha_minus_inv_beta*normal2* (face_surface);
992 robin_contribution= val * (marq ? tab_porosite_face(face) : 1) ;
993 tab_matrice(i, j) += robin_contribution;
994 tab_matrice(j ,i) += robin_contribution;
995 }
996 }
997 }
998 }
999 }
1000
1001 CDoubleArrView h_impose = static_cast<const ArrOfDouble&>(tab_h_impose).view_ro();
1002 CDoubleArrView derivee_flux_exterieur_imposee = static_cast<const ArrOfDouble&>(tab_derivee_flux_exterieur_imposee).view_ro();
1003 CDoubleArrView face_surfaces = domaine_VEF.face_surfaces().view_ro();
1004
1005 // Neumann: calcul des contributions aux bords
1006 auto neumann = KOKKOS_LAMBDA (const int face)
1007 {
1008 double h = h_impose(face);
1009 double dphi_dT = derivee_flux_exterieur_imposee(face);
1010 matrice.add(face,face, + (h + dphi_dT) * face_surfaces(face));
1011 };
1012
1013 Kokkos::parallel_for(start_gpu_timer(__KERNEL_NAME__), premiere_face_int, neumann);
1014 end_gpu_timer(__KERNEL_NAME__);
1015
1016 modifier_matrice_pour_periodique_apres_contribuer(tab_matrice,equation());
1017}
1018
1019void Op_Diff_VEF_Face::ajouter_contribution_multi_scalaire(const DoubleTab& tab_transporte, Matrice_Morse& tab_matrice) const
1020{
1021 modifier_matrice_pour_periodique_avant_contribuer(tab_matrice, equation());
1022
1023 // On remplit le tableau nu car l'assemblage d'une matrice
1024 // avec ajouter_contribution peut se faire avant le premier
1025 // pas de temps
1026 remplir_nu(nu_);
1027
1028 const Domaine_Cl_VEF& domaine_Cl_VEF = la_zcl_vef.valeur();
1029 const Domaine_VEF& domaine_VEF = le_dom_vef.valeur();
1030 const IntTab& tab_elem_faces = domaine_VEF.elem_faces();
1031 const IntTab& tab_face_voisins = domaine_VEF.face_voisins();
1032 const ArrOfInt& tab_est_face_bord = domaine_VEF.est_face_bord();
1033
1034 int nb_dim = tab_transporte.nb_dim();
1035 int nb_comp = (nb_dim == 2 ? tab_transporte.dimension(1) : 1);
1036
1037 DoubleTab tab_nu;
1038 int marq = phi_psi_diffuse(equation());
1039 const DoubleVect& porosite_elem = equation().milieu().porosite_elem();
1040
1041 // soit on a div(phi nu grad inco)
1042 // soit on a div(nu grad phi inco)
1043 // cela depend si on diffuse phi_psi ou psi
1044 modif_par_porosite_si_flag(nu_, tab_nu, !marq, porosite_elem);
1045 DoubleVect tab_porosite_eventuelle(equation().milieu().porosite_face());
1046 if (!marq)
1047 tab_porosite_eventuelle = 1;
1048
1049 int nb_faces_elem = domaine_VEF.domaine().nb_faces_elem();
1050 int nb_bords = domaine_VEF.nb_front_Cl();
1051 int nb_faces_tot = domaine_VEF.nb_faces_tot();
1052
1053 IntVect tab_face_associee(domaine_VEF.premiere_face_int());
1054 IntVect tab_fac2b_idx(domaine_VEF.nb_faces_tot());
1055
1056 // Recuperer les indices des faces de bord periodiques et les
1057 // faces associees dans des tableaux au prealable pour acces
1058 // structure dans le kernel ensuite
1059 for (int n_bord = 0; n_bord < nb_bords; n_bord++)
1060 {
1061 const Cond_lim& la_cl = domaine_Cl_VEF.les_conditions_limites(n_bord);
1062 const Front_VF& le_bord = ref_cast(Front_VF,la_cl->frontiere_dis());
1063 int num1 = le_bord.num_premiere_face();
1064
1065 if (sub_type(Periodique, la_cl.valeur()))
1066 {
1067 const Periodique& la_cl_perio = ref_cast(Periodique, la_cl.valeur());
1068
1069 int nb_faces = le_bord.nb_faces();
1070 int num2b = num1 + nb_faces / 2;
1071
1072 // on ne parcourt que la moitie des faces periodiques
1073 // on copiera a la fin le resultat dans la face
1074 // associee...
1075 ToDo_Kokkos("critical");
1076 for (int fac = num1; fac < num2b; fac++)
1077 {
1078 int fac_asso = la_cl_perio.face_associee(fac - num1) + num1;
1079
1080 tab_face_associee(fac) = fac_asso;
1081 tab_fac2b_idx(fac) = fac;
1082 tab_fac2b_idx(fac+nb_faces/2) = fac;
1083 }
1084 }
1085 }
1086
1087 CIntArrView est_face_bord = tab_est_face_bord.view_ro();
1088 CIntTabView face_voisins = tab_face_voisins.view_ro();
1089 CDoubleTabView face_normales = domaine_VEF.face_normales().view_ro();
1090 CDoubleArrView inverse_volumes = domaine_VEF.inverse_volumes().view_ro();
1091 CIntTabView elem_faces = tab_elem_faces.view_ro();
1092
1093 CDoubleTabView nu = tab_nu.view_ro();
1094 CDoubleArrView porosite_eventuelle = tab_porosite_eventuelle.view_ro();
1095
1096 CIntArrView fac2b_idx = tab_fac2b_idx.view_ro();
1097 CIntArrView face_associee = tab_face_associee.view_ro();
1098
1099 Matrice_Morse_View matrice;
1100 matrice.set(tab_matrice);
1101
1102 auto kern_elem_faces = KOKKOS_LAMBDA (const int fac)
1103 {
1104 int type_face = est_face_bord(fac);
1105 int fac2b = fac2b_idx(fac);
1106
1107 // Periodic boundary faces
1108 if (type_face == 2 && fac == fac2b)
1109 {
1110 int fac_asso = face_associee(fac);
1111
1112 int elem1 = face_voisins(fac, 0);
1113 int elem2 = face_voisins(fac, 1);
1114
1115 for (int i = 0; i < nb_faces_elem; i++)
1116 {
1117 int j = elem_faces(elem1, i);
1118 if (j > fac)
1119 {
1120 for (int nc = 0; nc < nb_comp; nc++)
1121 {
1122 double val = viscA(fac, j, elem1, nu(elem1, nc), face_voisins, face_normales, inverse_volumes);
1123
1124 int n0 = fac * nb_comp + nc;
1125 int j0 = j * nb_comp + nc;
1126
1127 matrice.atomic_add(n0, n0, + val * porosite_eventuelle(fac));
1128 matrice.atomic_add(n0, j0, - val * porosite_eventuelle(j));
1129 matrice.atomic_add(j0, n0, - val * porosite_eventuelle(fac));
1130 matrice.atomic_add(j0, j0, + val * porosite_eventuelle(j));
1131 }
1132 }
1133 if (elem2 != -1)
1134 {
1135 j = elem_faces(elem2, i);
1136 if (j > fac)
1137 {
1138 for (int nc = 0; nc < nb_comp; nc++)
1139 {
1140 double val = viscA(fac, j, elem2, nu(elem1, nc), face_voisins, face_normales, inverse_volumes);
1141
1142 int n0 = fac * nb_comp + nc;
1143 int j0 = j * nb_comp + nc;
1144 int n0perio = fac_asso * nb_comp + nc;
1145
1146 matrice.atomic_add(n0, n0, + val * porosite_eventuelle(fac));
1147 matrice.atomic_add(n0, j0, - val * porosite_eventuelle(j));
1148 matrice.atomic_add(j0, n0perio, - val * porosite_eventuelle(fac));
1149 matrice.atomic_add(j0, j0, + val * porosite_eventuelle(j));
1150 }
1151 }
1152 }
1153 }
1154 }
1155 // Faces de bord non periodiques
1156 else if (type_face == 1)
1157 {
1158 int elem1 = face_voisins(fac, 0);
1159
1160 for (int i = 0; i < nb_faces_elem; i++)
1161 {
1162 int j = elem_faces(elem1, i);
1163 if (j > fac)
1164 {
1165 for (int nc = 0; nc < nb_comp; nc++)
1166 {
1167 double val = viscA(fac, j, elem1, nu(elem1, nc), face_voisins, face_normales, inverse_volumes);
1168
1169 int n0 = fac * nb_comp + nc;
1170 int j0 = j * nb_comp + nc;
1171
1172 matrice.atomic_add(n0, n0, + val * porosite_eventuelle(fac));
1173 matrice.atomic_add(n0, j0, - val * porosite_eventuelle(j));
1174 matrice.atomic_add(j0, n0, - val * porosite_eventuelle(fac));
1175 matrice.atomic_add(j0, j0, + val * porosite_eventuelle(j));
1176 }
1177 }
1178 }
1179 }
1180 // Faces internes
1181 else if (type_face == 0)
1182 {
1183 int elem1 = face_voisins(fac, 0);
1184 int elem2 = face_voisins(fac, 1);
1185
1186 for (int i = 0; i < nb_faces_elem; i++)
1187 {
1188 int j = elem_faces(elem1, i);
1189 if (j > fac)
1190 {
1191 for (int nc = 0; nc < nb_comp; nc++)
1192 {
1193 double val = viscA(fac, j, elem1, nu(elem1, nc), face_voisins, face_normales, inverse_volumes);
1194
1195 int n0 = fac * nb_comp + nc;
1196 int j0 = j * nb_comp + nc;
1197
1198 matrice.atomic_add(n0, n0, + val * porosite_eventuelle(fac));
1199 matrice.atomic_add(n0, j0, - val * porosite_eventuelle(j));
1200 matrice.atomic_add(j0, n0, - val * porosite_eventuelle(fac));
1201 matrice.atomic_add(j0, j0, + val * porosite_eventuelle(j));
1202 }
1203 }
1204
1205 if (elem2 != -1)
1206 {
1207 j = elem_faces(elem2, i);
1208 if (j > fac)
1209 {
1210 for (int nc = 0; nc < nb_comp; nc++)
1211 {
1212 double val = viscA(fac, j, elem2, nu(elem2, nc), face_voisins, face_normales, inverse_volumes);
1213 int n0 = fac * nb_comp + nc;
1214 int j0 = j * nb_comp + nc;
1215
1216 matrice.atomic_add(n0, n0, + val * porosite_eventuelle(fac));
1217 matrice.atomic_add(n0, j0, - val * porosite_eventuelle(j));
1218 matrice.atomic_add(j0, n0, - val * porosite_eventuelle(fac));
1219 matrice.atomic_add(j0, j0, + val * porosite_eventuelle(j));
1220 }
1221 }
1222 }
1223 }
1224 }
1225 };
1226
1227 Kokkos::parallel_for(start_gpu_timer(__KERNEL_NAME__), nb_faces_tot, kern_elem_faces);
1228 end_gpu_timer(__KERNEL_NAME__);
1229
1230 for (int n_bord = 0; n_bord < nb_bords; n_bord++)
1231 {
1232 const Cond_lim& la_cl = domaine_Cl_VEF.les_conditions_limites(n_bord);
1233 const Front_VF& le_bord = ref_cast(Front_VF,la_cl->frontiere_dis());
1234
1235 if (sub_type(Echange_externe_impose, la_cl.valeur()))
1236 {
1237 const Echange_externe_impose& la_cl_paroi = ref_cast(Echange_externe_impose, la_cl.valeur());
1238 int ndeb = le_bord.num_premiere_face();
1239 int nfin = ndeb + le_bord.nb_faces();
1240 ToDo_Kokkos("critical");
1241 for (int face = ndeb; face < nfin; face++)
1242 for (int nc = 0; nc < nb_comp; nc++)
1243 {
1244 const int i = face * nb_comp + nc;
1245 tab_matrice(i, i) += la_cl_paroi.h_imp(face - ndeb, nc) * domaine_VEF.surface(face);
1246 }
1247 }
1248 if (sub_type(Echange_externe_radiatif, la_cl.valeur()))
1249 {
1250 throw;
1251 }
1252 }
1253
1254 modifier_matrice_pour_periodique_apres_contribuer(tab_matrice, equation());
1255}
1256
1257void Op_Diff_VEF_Face::contribue_au_second_membre(DoubleTab& resu ) const
1258{
1259 const Domaine_Cl_VEF& domaine_Cl_VEF = la_zcl_vef.valeur();
1260 const Domaine_VEF& domaine_VEF = le_dom_vef.valeur();
1261 int nb_comp = 1;
1262 int nb_dim = resu.nb_dim();
1263
1264 int nb_bords=domaine_VEF.nb_front_Cl();
1265
1266 if(nb_dim==2)
1267 nb_comp=resu.dimension(1);
1268
1269 // Partie imposee :
1270
1271 if (nb_dim == 1)
1272 {
1273 for (int n_bord=0; n_bord<nb_bords; n_bord++)
1274 {
1275 const Cond_lim& la_cl = domaine_Cl_VEF.les_conditions_limites(n_bord);
1276
1277 if (sub_type(Neumann_paroi,la_cl.valeur()))
1278 {
1279 const Neumann_paroi& la_cl_paroi = ref_cast(Neumann_paroi, la_cl.valeur());
1280 const Front_VF& le_bord = ref_cast(Front_VF,la_cl->frontiere_dis());
1281 int ndeb = le_bord.num_premiere_face();
1282 int nfin = ndeb + le_bord.nb_faces();
1283 ToDo_Kokkos("critical");
1284 for (int face=ndeb; face<nfin; face++)
1285 resu[face] += la_cl_paroi.flux_impose(face-ndeb)*domaine_VEF.surface(face);
1286 }
1287 else if (sub_type(Echange_externe_impose,la_cl.valeur()))
1288 {
1289 Cerr << "Non code pour Echange_externe_impose" << finl;
1290 assert(0);
1291 }
1292 else if (sub_type(Echange_externe_radiatif,la_cl.valeur()))
1293 {
1294 Cerr << "Non code pour Echange_externe_radiatif" << finl;
1295 assert(0);
1296 }
1297
1298 }
1299 }
1300 else
1301 {
1302 for (int n_bord=0; n_bord<nb_bords; n_bord++)
1303 {
1304 const Cond_lim& la_cl = domaine_Cl_VEF.les_conditions_limites(n_bord);
1305
1306 if (sub_type(Neumann_paroi,la_cl.valeur()))
1307 {
1308 const Neumann_paroi& la_cl_paroi = ref_cast(Neumann_paroi, la_cl.valeur());
1309 const Front_VF& le_bord = ref_cast(Front_VF,la_cl->frontiere_dis());
1310 int ndeb = le_bord.num_premiere_face();
1311 int nfin = ndeb + le_bord.nb_faces();
1312 ToDo_Kokkos("critical");
1313 for (int face=ndeb; face<nfin; face++)
1314 for (int comp=0; comp<nb_comp; comp++)
1315 resu(face,comp) += la_cl_paroi.flux_impose(face-ndeb,comp)*domaine_VEF.surface(face);
1316 }
1317 }
1318 }
1319}
1320
1321void Op_Diff_VEF_Face::verifier() const
1322{
1323 static int testee=0;
1324 if(testee)
1325 return;
1326 testee=1;
1327 const Domaine_VEF& domaine_VEF = le_dom_vef.valeur();
1328 // const Domaine_Cl_VEF& domaine_Cl_VEF = la_zcl_vef.valeur();
1329 // const Conds_lim& les_cl = domaine_Cl_VEF.les_conditions_limites();
1330 const DoubleVect& volumes_entrelaces = domaine_VEF.volumes_entrelaces();
1331
1332 const DoubleTab& xv=domaine_VEF.xv();
1333 DoubleTab vit(equation().inconnue().valeurs());
1334 DoubleTab resu(vit);
1335 int i, comp;
1336 if(dimension==2)
1337 {
1338 const int nbf = vit.dimension(0);
1339 Cerr << " Verification de delta(x,0) " << finl;
1340 for(i=0; i<nbf; i++)
1341 {
1342 vit(i,0)=xv(i,0);
1343 vit(i,1)=0;
1344 }
1345 calculer(vit, resu);
1346 for(i=0; i<nbf; i++)
1347 for(comp=0; comp<dimension; comp++)
1348 resu(i,comp)/=(volumes_entrelaces(i));
1349 for(i=0; i<nbf; i++)
1350 {
1351 if(std::fabs(resu(i,0))>1.e-10)
1352 {
1353 Cerr << " delta(x,0) ("<<i<<") = "
1354 << resu(i,0);
1355 Cerr << finl;
1356 }
1357 }
1358 Cerr << " Verification de delta(y(1-y),0) " << finl;
1359 for(i=0; i<nbf; i++)
1360 {
1361 vit(i,0)=xv(i,1)*(1-xv(i,1));
1362 vit(i,1)=0;
1363 }
1364 calculer(vit, resu);
1365 for(i=0; i<nbf; i++)
1366 for(comp=0; comp<dimension; comp++)
1367 resu(i,comp)/=(volumes_entrelaces(i));
1368 for(i=0; i<nbf; i++)
1369 {
1370 if(std::fabs(2-resu(i,0))>1.e-10)
1371 {
1372 Cerr << " delta(y(1-y),0) ("<<i<<") = "
1373 << resu(i,0);
1374 Cerr << finl;
1375 }
1376 }
1377 }
1378}
1379
1380
1381
Champ_Inc_base::valeurs
DoubleTab & valeurs() override
Renvoie le tableau des valeurs du champ au temps courant.
Definition Champ_Inc_base.h:77
Champ_P1NC::calcul_gradient
static DoubleTab & calcul_gradient(const DoubleTab &, DoubleTab &, const Domaine_Cl_VEF &)
Definition Champ_P1NC.cpp:915
Champ_base
classe Champ_base Cette classe est la base de la hierarchie des champs.
Definition Champ_base.h:43
Cond_lim
classe Cond_lim Classe generique servant a representer n'importe quelle classe
Definition Cond_lim.h:31
Domaine_32_64::creer_tableau_elements
virtual void creer_tableau_elements(Array_base &, RESIZE_OPTIONS opt=RESIZE_OPTIONS::COPY_INIT) const
creation d'un tableau parallele de valeurs aux elements.
Definition Domaine.cpp:851
Domaine_32_64::nb_faces_elem
int nb_faces_elem(int=0) const
Renvoie le nombre de face de type i des elements geometriques constituants le domaine.
Definition Domaine.h:484
Domaine_Cl_VEF
Definition Domaine_Cl_VEF.h:40
Domaine_Cl_dis_base::les_conditions_limites
const Cond_lim & les_conditions_limites(int) const
Renvoie la i-ieme condition aux limites.
Definition Domaine_Cl_dis_base.cpp:386
Domaine_VEF
class Domaine_VEF
Definition Domaine_VEF.h:54
Domaine_VF::face_surfaces
virtual const DoubleVect & face_surfaces() const
Definition Domaine_VF.h:51
Domaine_VF::nb_faces
int nb_faces() const
renvoie le nombre global de faces.
Definition Domaine_VF.h:471
Domaine_VF::volumes_entrelaces
DoubleVect & volumes_entrelaces()
Definition Domaine_VF.h:99
Domaine_VF::nb_faces_tot
int nb_faces_tot() const
renvoie le nombre total de faces.
Definition Domaine_VF.h:481
Domaine_VF::face_normales
virtual double face_normales(int face, int comp) const
Definition Domaine_VF.h:47
Domaine_VF::xv
double xv(int num_face, int k) const
Definition Domaine_VF.h:76
Domaine_VF::surface
virtual double surface(int i) const
Definition Domaine_VF.h:53
Domaine_VF::normalized_boundaries_outward_vector
DoubleTab normalized_boundaries_outward_vector(int global_face_number, double scale_factor) const
Compute the normalized boundary outward vector associated to the face global_face_number and eventual...
Definition Domaine_VF.cpp:717
Domaine_VF::elem_faces
int elem_faces(int i, int j) const
renvoie le numero de le ieme face de la maille num_elem la facon dont ces faces sont numerotees est
Definition Domaine_VF.h:543
Domaine_VF::est_face_bord
ArrOfInt & est_face_bord()
Definition Domaine_VF.h:85
Domaine_VF::premiere_face_int
int premiere_face_int() const
une face est interne ssi elle separe deux elements.
Definition Domaine_VF.h:463
Domaine_VF::face_voisins
int face_voisins(int num_face, int i) const
renvoie l'element voisin de numface dans la direction i.
Definition Domaine_VF.h:418
Domaine_VF::nb_faces_bord
int nb_faces_bord() const
renvoie le nombre de faces sur lesquelles sont appliquees les conditions limites :
Definition Domaine_VF.h:513
Domaine_VF::inverse_volumes
double inverse_volumes(int i) const
Definition Domaine_VF.h:114
Domaine_dis_base::nb_front_Cl
int nb_front_Cl() const
Definition Domaine_dis_base.h:53
Domaine_dis_base::domaine
const Domaine & domaine() const
Definition Domaine_dis_base.h:46
Echange_couplage_thermique
: class Echange_couplage_thermique
Definition Echange_couplage_thermique.h:26
Echange_couplage_thermique::flux_exterieur_impose
double flux_exterieur_impose(int i) const override
Definition Echange_couplage_thermique.h:40
Echange_externe_impose
Classe Echange_externe_impose: Cette classe represente le cas particulier de la classe.
Definition Echange_externe_impose.h:46
Echange_externe_radiatif
classe Echange_externe_radiatif: Combines radiative (sigma * eps * (T^4 - T_ext^4)) and convective (h...
Definition Echange_externe_radiatif.h:31
Echange_global_impose::derivee_flux_exterieur_imposee
virtual double derivee_flux_exterieur_imposee(int i) const
Definition Echange_global_impose.cpp:244
Echange_impose_base::h_imp
virtual double h_imp(int num) const
Renvoie la valeur du coefficient d'echange de chaleur impose sur la i-eme composante.
Definition Echange_impose_base.cpp:112
Echange_impose_base::has_emissivite
bool has_emissivite() const
Definition Echange_impose_base.h:45
Echange_impose_base::T_ext
virtual double T_ext(int num) const
Renvoie la valeur de la temperature imposee sur la i-eme composante du champ de frontiere.
Definition Echange_impose_base.cpp:76
Echange_impose_base::tab_h_imp
const DoubleTab & tab_h_imp(double temps=DMAXFLOAT) const
Definition Echange_impose_base.cpp:284
Echange_impose_base::emissivite
double emissivite(int num) const
Renvoie la valeur de l'emissivite impose sur la i-eme composante.
Definition Echange_impose_base.cpp:150
Echange_impose_base::tab_T_ext
const DoubleTab & tab_T_ext(double temps=DMAXFLOAT) const
Definition Echange_impose_base.cpp:258
Echange_impose_base::tab_emissivite
const DoubleTab & tab_emissivite(double temps=DMAXFLOAT) const
Definition Echange_impose_base.cpp:310
Entree
Class defining operators and methods for all reading operation in an input flow (file,...
Definition Entree.h:42
Equation_base::milieu
virtual const Milieu_base & milieu() const =0
Equation_base::inconnue
virtual const Champ_Inc_base & inconnue() const =0
Field_base::nature_du_champ
virtual Nature_du_champ nature_du_champ() const
Definition Field_base.h:77
Front_VF
class Front_VF
Definition Front_VF.h:36
Front_VF::nb_faces
int nb_faces() const
Definition Front_VF.h:53
Front_VF::num_premiere_face
int num_premiere_face() const
Definition Front_VF.h:63
Front_VF::nb_faces_tot
int nb_faces_tot() const
Definition Front_VF.h:58
Front_VF::num_face
int num_face(const int) const
Definition Front_VF.h:68
Matrice_Morse
Classe Matrice_Morse Represente une matrice M (creuse), non necessairement carree.
Definition Matrice_Morse.h:50
Milieu_base::porosite_elem
DoubleVect & porosite_elem()
Definition Milieu_base.h:58
Milieu_base::porosite_face
DoubleVect & porosite_face()
Definition Milieu_base.h:62
MorEqn::equation
const Equation_base & equation() const
Renvoie la reference sur l'equation pointe par MorEqn::mon_equation.
Definition MorEqn.h:62
Neumann_homogene
Classe Neumann_homogene Cette classe est la classe de base de la hierarchie des conditions aux limite...
Definition Neumann_homogene.h:31
Neumann_paroi
Classe Neumann_paroi Cette condition limite correspond a un flux impose pour l'equation de.
Definition Neumann_paroi.h:29
Neumann_sortie_libre
classe Neumann_sortie_libre Cette classe represente une frontiere ouverte sans vitesse imposee
Definition Neumann_sortie_libre.h:34
Neumann::flux_impose
virtual double flux_impose(int i) const
Renvoie la valeur du flux impose sur la i-eme composante du champ representant le flux a la frontiere...
Definition Neumann.cpp:35
Objet_U::dimension
static int dimension
Definition Objet_U.h:99
Objet_U::que_suis_je
const Nom & que_suis_je() const
renvoie la chaine identifiant la classe.
Definition Objet_U.cpp:104
Objet_U::readOn
virtual Entree & readOn(Entree &)
Lecture d'un Objet_U sur un flot d'entree Methode a surcharger.
Definition Objet_U.cpp:293
Objet_U::printOn
virtual Sortie & printOn(Sortie &) const
Ecriture de l'objet sur un flot de sortie Methode a surcharger.
Definition Objet_U.cpp:282
Op_Diff_VEF_Face
class Op_Diff_VEF_Face Cette classe represente l'operateur de diffusion
Definition Op_Diff_VEF_Face.h:32
Op_Diff_VEF_Face::diffusivite
const Champ_base & diffusivite() const override
Definition Op_Diff_VEF_Face.cpp:66
Op_Diff_VEF_Face::Op_Diff_VEF_Face
Op_Diff_VEF_Face()
Definition Op_Diff_VEF_Face.cpp:38
Op_Diff_VEF_Face::ajouter
DoubleTab & ajouter(const DoubleTab &, DoubleTab &) const override
Definition Op_Diff_VEF_Face.cpp:679
Op_Diff_VEF_Face::contribue_au_second_membre
void contribue_au_second_membre(DoubleTab &) const
Definition Op_Diff_VEF_Face.cpp:1257
Op_Diff_VEF_Face::ajouter_cas_scalaire
void ajouter_cas_scalaire(const DoubleTab &inconnue, DoubleTab &resu, DoubleTab &flux_bords, DoubleTab &nu, const Domaine_Cl_VEF &domaine_Cl_VEF, const Domaine_VEF &domaine_VEF) const
Definition Op_Diff_VEF_Face.cpp:76
Op_Diff_VEF_Face::calculer
DoubleTab & calculer(const DoubleTab &, DoubleTab &) const override
Definition Op_Diff_VEF_Face.cpp:718
Op_Diff_VEF_Face::completer
void completer() override
Associe l'operateur au domaine_dis, le domaine_Cl_dis, et a l'inconnue de son equation.
Definition Op_Diff_VEF_Face.cpp:61
Op_Diff_VEF_Face::ajouter_cas_vectoriel
void ajouter_cas_vectoriel(const DoubleTab &inconnue, DoubleTab &resu, DoubleTab &flux_bords, DoubleTab &nu, const Domaine_Cl_VEF &domaine_Cl_VEF, const Domaine_VEF &domaine_VEF, int nb_comp) const
Definition Op_Diff_VEF_Face.cpp:291
Op_Diff_VEF_Face::grad_
DoubleTab grad_
Definition Op_Diff_VEF_Face.h:78
Op_Diff_VEF_Face::verifier
void verifier() const
Definition Op_Diff_VEF_Face.cpp:1321
Op_Diff_VEF_Face::ajouter_cas_multi_scalaire
void ajouter_cas_multi_scalaire(const DoubleTab &inconnue, DoubleTab &resu, DoubleTab &flux_bords, DoubleTab &nu, const Domaine_Cl_VEF &domaine_Cl_VEF, const Domaine_VEF &domaine_VEF, int nb_comp) const
Definition Op_Diff_VEF_Face.cpp:461
Op_Diff_VEF_Face::ajouter_contribution_multi_scalaire
void ajouter_contribution_multi_scalaire(const DoubleTab &, Matrice_Morse &) const
Definition Op_Diff_VEF_Face.cpp:1019
Op_Diff_VEF_Face::associer_diffusivite
void associer_diffusivite(const Champ_base &) override
associe le champ de diffusivite
Definition Op_Diff_VEF_Face.cpp:56
Op_Diff_VEF_Face::ajouter_contribution
void ajouter_contribution(const DoubleTab &, Matrice_Morse &) const
Definition Op_Diff_VEF_Face.cpp:724
Op_Diff_VEF_base
class Op_Diff_VEF_base
Definition Op_Diff_VEF_base.h:39
Op_Diff_VEF_base::phi_psi_diffuse
int phi_psi_diffuse(const Equation_base &eq) const
definit si on calcule div(phi nu grad Psi) ou div(nu grap Phi psi)
Definition Op_Diff_VEF_base.cpp:44
Op_Diff_VEF_base::nu_
DoubleTab nu_
Definition Op_Diff_VEF_base.h:66
Op_Diff_VEF_base::remplir_nu
virtual void remplir_nu(DoubleTab &) const
Definition Op_Diff_VEF_base.cpp:251
Op_Diff_VEF_base::viscA
double viscA(int face_i, int face_j, int num_elem, const _TYPE_ &diffu) const
Definition Op_Diff_VEF_base.h:92
Op_VEF_Face::modifier_matrice_pour_periodique_apres_contribuer
void modifier_matrice_pour_periodique_apres_contribuer(Matrice_Morse &matrice, const Equation_base &) const
Somme les 2 lignes des faces periodiques associees permet de calculer dans le code sans se poser de q...
Definition Op_VEF_Face.cpp:738
Op_VEF_Face::modifier_matrice_pour_periodique_avant_contribuer
void modifier_matrice_pour_periodique_avant_contribuer(Matrice_Morse &matrice, const Equation_base &) const
divise les coefficients sur les ligne des faces periodiques par 2 en prevision de l'application modif...
Definition Op_VEF_Face.cpp:725
Op_VEF_Face::modifier_flux
void modifier_flux(const Operateur_base &) const
Definition Op_VEF_Face.cpp:338
Operateur_base::flux_bords_
DoubleTab flux_bords_
Definition Operateur_base.h:140
Operateur_base::completer
virtual void completer()
Associe l'operateur au domaine_dis, le domaine_Cl_dis, et a l'inconnue de son equation.
Definition Operateur_base.cpp:89
Operateur_base::flux_bords
DoubleTab & flux_bords()
Definition Operateur_base.h:107
Periodique
classe Periodique Cette classe represente une condition aux limites periodique.
Definition Periodique.h:31
Periodique::face_associee
int face_associee(int i) const
Definition Periodique.h:35
Robin_VEF
Class Robin_VEF for Robin boundary conditions.
Definition Robin_VEF.h:44
Robin_VEF::get_beta_cl
double get_beta_cl() const
Definition Robin_VEF.h:53
Robin_VEF::flux_tangentiel_imp
double flux_tangentiel_imp(int i, int j) const
Returns the value of the imposed tangential flux on the (i,j)-th component of the field representing ...
Definition Robin_VEF.cpp:111
Robin_VEF::flux_normal_imp
double flux_normal_imp(int i) const
Returns the value of the imposed normal flux on the i-th component of the field representing the flux...
Definition Robin_VEF.cpp:103
Robin_VEF::get_alpha_cl
double get_alpha_cl() const
Definition Robin_VEF.h:52
Sortie
Classe de base des flux de sortie.
Definition Sortie.h:52
Support_Champ_Masse_Volumique::declare_support_masse_volumique
virtual void declare_support_masse_volumique(int ok)
Le constructeur d'une classe derivee qui se sert de la masse volumique doit appeler cette fonction av...
Definition Support_Champ_Masse_Volumique.cpp:46
Symetrie
classe Symetrie Sur les faces de symetrie on a les proprietes suivantes:
Definition Symetrie.h:37
TRUSTTab::view_wo
std::enable_if_t< is_default_exec_space< EXEC_SPACE >, View< _TYPE_, _SHAPE_ > > view_wo()
Definition TRUSTTab.h:276
TRUSTTab::nb_dim
int nb_dim() const
Definition TRUSTTab.h:199
TRUSTTab::view_ro
std::enable_if_t< is_default_exec_space< EXEC_SPACE >, ConstView< _TYPE_, _SHAPE_ > > view_ro() const
Definition TRUSTTab.h:261
TRUSTTab::view_rw
std::enable_if_t< is_default_exec_space< EXEC_SPACE >, View< _TYPE_, _SHAPE_ > > view_rw()
Definition TRUSTTab.h:291
TRUSTTab::dimension
_SIZE_ dimension(int d) const
Definition TRUSTTab.tpp:133