| (1) assembly->set_expression(elements_expression |
| (2) ( |
| (3) AllElementsT(), |
| (4) group |
| (5) ( |
| (6) _A = _0, _T = _0, |
| (7) compute_tau(u, nu_eff, u_ref, lit(tau_ps), lit(tau_su), lit(tau_bulk)), |
| (8) element_quadrature |
| (9) ( |
| (10) _A(p, u[_i]) += transpose(N(p) + tau_ps*u_adv*nabla(p)*0.5) * nabla(u)[_i] |
| (11) + tau_ps * transpose(nabla(p)[_i]) * u_adv*nabla(u), |
| (12) _A(p, p) += tau_ps * transpose(nabla(p)) * nabla(p) / rho, |
| (13) _A(u[_i], u[_j]) += transpose((tau_bulk + 1/3*nu_eff)*nabla(u)[_i] |
| (14) + 0.5*u_adv[_i]*(N(u) + tau_su*u_adv*nabla(u))) * nabla(u)[_j], |
| (15) _A(u[_i], u[_i]) += nu_eff * transpose(nabla(u)) * nabla(u) |
| (16) + transpose(N(u) + tau_su*u_adv*nabla(u)) * u_adv*nabla(u), |
| (17) _A(u[_i], p) += transpose(N(u) + tau_su*u_adv*nabla(u)) * nabla(p)[_i] / rho, |
| (18) _T(p, u[_i]) += tau_ps * transpose(nabla(p)[_i]) * N(u), |
| (19) _T(u[_i], u[_i]) += transpose(N(u) + tau_su*u_adv*nabla(u)) * N(u) |
| (20) ), |
| (21) system_rhs += −_A * _x, |
| (22) _A(p) = _A(p) / theta, |
| (23) system_matrix += invdt() * _T + theta * _A |
| (24) ))); |
