370 {
372
373 const double vol = getMeasure();
374 auto t_dot_u = getFTensor1FromMat<U_FIELD_DIM>(*
dotUPtr);
375 auto t_u = getFTensor1FromMat<U_FIELD_DIM>(*
uPtr);
376 auto t_p = getFTensor0FromVec(*
pPtr);
377 auto t_grad_u = getFTensor2FromMat<U_FIELD_DIM, SPACE_DIM>(*
gradUPtr);
378 auto t_h = getFTensor0FromVec(*
hPtr);
379 auto t_grad_h = getFTensor1FromMat<SPACE_DIM>(*
gradHPtr);
380 auto t_g = getFTensor0FromVec(*
gPtr);
381 auto t_coords = getFTensor1CoordsAtGaussPts();
382
383 auto t_base = data.getFTensor0N();
384 auto t_diff_base = data.getFTensor1DiffN<
SPACE_DIM>();
385
386 auto t_w = getFTensor0IntegrationWeight();
387
396
399
400 for (int gg = 0; gg != nbIntegrationPts; gg++) {
401
402 const double r = t_coords(0);
404
407
409
410 t_inertia_force(
i) = (
rho * alpha) * (t_dot_u(
i));
411
413 t_phase_force(
i) = -alpha *
kappa * t_g * t_grad_h(
i);
414 t_convection(
i) = (
rho * alpha) * (t_u(
j) * t_grad_u(
i,
j));
415
417 alpha * (t_D(
i,
j,
k,
l) * t_grad_u(
k,
l) +
t_kd(
i,
j) * t_p);
418
419 auto t_nf = getFTensor1FromArray<U_FIELD_DIM, U_FIELD_DIM>(locF);
420
421 t_forces(
i) = t_inertia_force(
i) + t_buoyancy(
i) + t_gravity(
i) +
422 t_convection(
i) + t_phase_force(
i);
423
424 int bb = 0;
426
427 t_nf(
i) += t_base * t_forces(
i);
428 t_nf(
i) += t_diff_base(
j) * t_stress(
i,
j);
429
431 t_nf(0) += (t_base * (alpha / t_coords(0))) * (2 *
mu * t_u(0) + t_p);
432 }
433
434 ++t_base;
435 ++t_diff_base;
436 ++t_nf;
437 }
438
439 for (; bb < nbRowBaseFunctions; ++bb) {
440 ++t_diff_base;
441 ++t_base;
442 }
443
444 ++t_dot_u;
445 ++t_u;
446 ++t_grad_u;
447 ++t_h;
448 ++t_grad_h;
449 ++t_g;
450 ++t_p;
451
452 ++t_w;
453 ++t_coords;
454 }
455
457 }
Kronecker Delta class symmetric.
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
constexpr int U_FIELD_DIM
FTensor::Index< 'i', SPACE_DIM > i
constexpr CoordinateTypes coord_type
FTensor::Index< 'l', 3 > l
FTensor::Index< 'j', 3 > j
FTensor::Index< 'k', 3 > k