A model for cardiac conduction. More...
#include <MMFDiffusion.h>
Inheritance diagram for Nektar::MMFDiffusion:
Static Public Member Functions | |
| static SolverUtils::EquationSystemSharedPtr | create (const LibUtilities::SessionReaderSharedPtr &pSession, const SpatialDomains::MeshGraphSharedPtr &pGraph) |
| Creates an instance of this class. More... | |
Static Public Attributes | |
| static std::string | className |
| Name of class. More... | |
 Static Public Attributes inherited from Nektar::SolverUtils::UnsteadySystem | |
| static std::string | cmdSetStartTime |
| static std::string | cmdSetStartChkNum |
Protected Attributes | |
| TestType | m_TestType |
| InitWaveType | m_InitWaveType |
| NekDouble | m_InitPtx |
| NekDouble | m_InitPty |
| NekDouble | m_InitPtz |
| Protected Attributes inherited from Nektar::SolverUtils::MMFSystem | |
| NekDouble | m_alpha |
| NekDouble | m_Incfreq |
| int | m_SmoothFactor |
| NekDouble | m_SFinit |
| Array< OneD, Array< OneD, NekDouble > > | m_movingframes |
| Array< OneD, Array< OneD, NekDouble > > | m_surfaceNormal |
| Array< OneD, Array< OneD, NekDouble > > | m_ncdotMFFwd |
| Array< OneD, Array< OneD, NekDouble > > | m_ncdotMFBwd |
| Array< OneD, Array< OneD, NekDouble > > | m_nperpcdotMFFwd |
| Array< OneD, Array< OneD, NekDouble > > | m_nperpcdotMFBwd |
| Array< OneD, Array< OneD, NekDouble > > | m_DivMF |
| Array< OneD, Array< OneD, Array< OneD, NekDouble > > > | m_CurlMF |
| Array< OneD, Array< OneD, Array< OneD, NekDouble > > > | m_MFtraceFwd |
| Array< OneD, Array< OneD, Array< OneD, NekDouble > > > | m_MFtraceBwd |
| Array< OneD, Array< OneD, Array< OneD, NekDouble > > > | m_ntimesMFFwd |
| Array< OneD, Array< OneD, Array< OneD, NekDouble > > > | m_ntimesMFBwd |
| Array< OneD, Array< OneD, Array< OneD, NekDouble > > > | m_ntimes_ntimesMFFwd |
| Array< OneD, Array< OneD, Array< OneD, NekDouble > > > | m_ntimes_ntimesMFBwd |
| Array< OneD, Array< OneD, NekDouble > > | m_ZimFwd |
| Array< OneD, Array< OneD, NekDouble > > | m_ZimBwd |
| Array< OneD, Array< OneD, NekDouble > > | m_YimFwd |
| Array< OneD, Array< OneD, NekDouble > > | m_YimBwd |
| Array< OneD, Array< OneD, NekDouble > > | m_epsvec |
| Array< OneD, Array< OneD, NekDouble > > | m_muvec |
| Array< OneD, Array< OneD, NekDouble > > | m_negepsvecminus1 |
| Array< OneD, Array< OneD, NekDouble > > | m_negmuvecminus1 |
| Array< OneD, Array< OneD, Array< OneD, Array< OneD, NekDouble > > > > | m_dedxi_cdot_e |
| SpatialDomains::GeomMMF | m_MMFdir |
| Protected Attributes inherited from Nektar::SolverUtils::UnsteadySystem | |
| LibUtilities::TimeIntegrationSchemeSharedPtr | m_intScheme |
| Wrapper to the time integration scheme. More... | |
| LibUtilities::TimeIntegrationSchemeOperators | m_ode |
| The time integration scheme operators to use. More... | |
| Array< OneD, Array< OneD, NekDouble > > | m_previousSolution |
| Storage for previous solution for steady-state check. More... | |
| std::vector< int > | m_intVariables |
| NekDouble | m_cflSafetyFactor |
| CFL safety factor (comprise between 0 to 1). More... | |
| NekDouble | m_CFLGrowth |
| CFL growth rate. More... | |
| NekDouble | m_CFLEnd |
| Maximun cfl in cfl growth. More... | |
| int | m_abortSteps |
| Number of steps between checks for abort conditions. More... | |
| bool | m_explicitDiffusion |
| Indicates if explicit or implicit treatment of diffusion is used. More... | |
| bool | m_explicitAdvection |
| Indicates if explicit or implicit treatment of advection is used. More... | |
| bool | m_explicitReaction |
| Indicates if explicit or implicit treatment of reaction is used. More... | |
| int | m_steadyStateSteps |
| Check for steady state at step interval. More... | |
| NekDouble | m_steadyStateTol |
| Tolerance to which steady state should be evaluated at. More... | |
| int | m_filtersInfosteps |
| Number of time steps between outputting filters information. More... | |
| std::vector< std::pair< std::string, FilterSharedPtr > > | m_filters |
| bool | m_homoInitialFwd |
| Flag to determine if simulation should start in homogeneous forward transformed state. More... | |
| std::ofstream | m_errFile |
| NekDouble | m_epsilon |
| Diffusion coefficient. More... | |
| Protected Attributes inherited from Nektar::SolverUtils::EquationSystem | |
| LibUtilities::CommSharedPtr | m_comm |
| Communicator. More... | |
| bool | m_verbose |
| LibUtilities::SessionReaderSharedPtr | m_session |
| The session reader. More... | |
| std::map< std::string, SolverUtils::SessionFunctionSharedPtr > | m_sessionFunctions |
| Map of known SessionFunctions. More... | |
| LibUtilities::FieldIOSharedPtr | m_fld |
| Field input/output. More... | |
| Array< OneD, MultiRegions::ExpListSharedPtr > | m_fields |
| Array holding all dependent variables. More... | |
| SpatialDomains::BoundaryConditionsSharedPtr | m_boundaryConditions |
| Pointer to boundary conditions object. More... | |
| SpatialDomains::MeshGraphSharedPtr | m_graph |
| Pointer to graph defining mesh. More... | |
| std::string | m_sessionName |
| Name of the session. More... | |
| NekDouble | m_time |
| Current time of simulation. More... | |
| int | m_initialStep |
| Number of the step where the simulation should begin. More... | |
| NekDouble | m_fintime |
| Finish time of the simulation. More... | |
| NekDouble | m_timestep |
| Time step size. More... | |
| NekDouble | m_lambda |
| Lambda constant in real system if one required. More... | |
| NekDouble | m_checktime |
| Time between checkpoints. More... | |
| NekDouble | m_lastCheckTime |
| NekDouble | m_TimeIncrementFactor |
| int | m_nchk |
| Number of checkpoints written so far. More... | |
| int | m_steps |
| Number of steps to take. More... | |
| int | m_checksteps |
| Number of steps between checkpoints. More... | |
| int | m_infosteps |
| Number of time steps between outputting status information. More... | |
| int | m_iterPIT = 0 |
| Number of parallel-in-time time iteration. More... | |
| int | m_windowPIT = 0 |
| Index of windows for parallel-in-time time iteration. More... | |
| int | m_spacedim |
| Spatial dimension (>= expansion dim). More... | |
| int | m_expdim |
| Expansion dimension. More... | |
| bool | m_singleMode |
| Flag to determine if single homogeneous mode is used. More... | |
| bool | m_halfMode |
| Flag to determine if half homogeneous mode is used. More... | |
| bool | m_multipleModes |
| Flag to determine if use multiple homogenenous modes are used. More... | |
| bool | m_useFFT |
| Flag to determine if FFT is used for homogeneous transform. More... | |
| bool | m_homogen_dealiasing |
| Flag to determine if dealiasing is used for homogeneous simulations. More... | |
| bool | m_specHP_dealiasing |
| Flag to determine if dealisising is usde for the Spectral/hp element discretisation. More... | |
| enum MultiRegions::ProjectionType | m_projectionType |
| Type of projection; e.g continuous or discontinuous. More... | |
| Array< OneD, Array< OneD, NekDouble > > | m_traceNormals |
| Array holding trace normals for DG simulations in the forwards direction. More... | |
| Array< OneD, bool > | m_checkIfSystemSingular |
| Flag to indicate if the fields should be checked for singularity. More... | |
| LibUtilities::FieldMetaDataMap | m_fieldMetaDataMap |
| Map to identify relevant solver info to dump in output fields. More... | |
| Array< OneD, NekDouble > | m_movingFrameData |
| Moving reference frame status in the inertial frame X, Y, Z, Theta_x, Theta_y, Theta_z, U, V, W, Omega_x, Omega_y, Omega_z, A_x, A_y, A_z, DOmega_x, DOmega_y, DOmega_z, pivot_x, pivot_y, pivot_z. More... | |
| std::vector< std::string > | m_strFrameData |
| variable name in m_movingFrameData More... | |
| int | m_NumQuadPointsError |
| Number of Quadrature points used to work out the error. More... | |
| enum HomogeneousType | m_HomogeneousType |
| NekDouble | m_LhomX |
| physical length in X direction (if homogeneous) More... | |
| NekDouble | m_LhomY |
| physical length in Y direction (if homogeneous) More... | |
| NekDouble | m_LhomZ |
| physical length in Z direction (if homogeneous) More... | |
| int | m_npointsX |
| number of points in X direction (if homogeneous) More... | |
| int | m_npointsY |
| number of points in Y direction (if homogeneous) More... | |
| int | m_npointsZ |
| number of points in Z direction (if homogeneous) More... | |
| int | m_HomoDirec |
| number of homogenous directions More... | |
| Protected Attributes inherited from Nektar::SolverUtils::ALEHelper | |
| Array< OneD, MultiRegions::ExpListSharedPtr > | m_fieldsALE |
| Array< OneD, Array< OneD, NekDouble > > | m_gridVelocity |
| Array< OneD, Array< OneD, NekDouble > > | m_gridVelocityTrace |
| std::vector< ALEBaseShPtr > | m_ALEs |
| bool | m_ALESolver = false |
| bool | m_ImplicitALESolver = false |
| NekDouble | m_prevStageTime = 0.0 |
| int | m_spaceDim |
Private Attributes | |
| StdRegions::VarCoeffMap | m_varcoeff |
| Variable diffusivity. More... | |
| Array< OneD, NekDouble > | m_epsilon |
| Array< OneD, NekDouble > | m_epsu |
Friends | |
| class | MemoryManager< MMFDiffusion > |
Detailed Description
A model for cardiac conduction.
Definition at line 79 of file MMFDiffusion.h.
Constructor & Destructor Documentation
◆ MMFDiffusion()
|
protected |
Constructor.
Definition at line 58 of file MMFDiffusion.cpp.
61{
62}
SOLVER_UTILS_EXPORT MMFSystem(const LibUtilities::SessionReaderSharedPtr &pSession, const SpatialDomains::MeshGraphSharedPtr &pGraph)
Definition: MMFSystem.cpp:41
SOLVER_UTILS_EXPORT UnsteadySystem(const LibUtilities::SessionReaderSharedPtr &pSession, const SpatialDomains::MeshGraphSharedPtr &pGraph)
Initialises UnsteadySystem class members.
Definition: UnsteadySystem.cpp:74
◆ ~MMFDiffusion()
|
overrideprotecteddefault |
Member Function Documentation
◆ create()
|
inlinestatic |
Creates an instance of this class.
Definition at line 85 of file MMFDiffusion.h.
88 {
90 MemoryManager<MMFDiffusion>::AllocateSharedPtr(pSession, pGraph);
91 p->InitObject();
93 }
static std::shared_ptr< DataType > AllocateSharedPtr(const Args &...args)
Allocate a shared pointer from the memory pool.
Definition: NekMemoryManager.hpp:166
std::shared_ptr< EquationSystem > EquationSystemSharedPtr
A shared pointer to an EquationSystem object.
Definition: EquationSystem.h:65
References Nektar::MemoryManager< DataType >::AllocateSharedPtr(), and CellMLToNektar.cellml_metadata::p.
◆ DoImplicitSolve()
|
protected |
Solve for the diffusion term.
OdeRhs
- Parameters
-
inarray Input array. outarray Output array. time Current simulation time. lambda Timestep.
Definition at line 185 of file MMFDiffusion.cpp.
189{
190 int nvariables = inarray.size();
192
195
196 Array<OneD, Array<OneD, NekDouble>> F(nvariables);
198 F[0] = Array<OneD, NekDouble>(nq * nvariables);
199
200 for (int n = 1; n < nvariables; ++n)
201 {
202 F[n] = F[n - 1] + nq;
203 }
204
205 // We solve ( \nabla^2 - HHlambda ) Y[i] = rhs [i]
206 // inarray = input: \hat{rhs} -> output: \hat{Y}
207 // outarray = output: nabla^2 \hat{Y}
208 // where \hat = modal coeffs
209 SetBoundaryConditions(time);
210
211 for (int i = 0; i < nvariables; ++i)
212 {
214
215 // Multiply 1.0/timestep
217 F[i], 1);
218
219 // Solve a system of equations with Helmholtz solver and transform
220 // back into physical space.
222 m_varcoeff);
223
225 }
226}
Definition: BasicUtils/SharedArray.hpp:57
Array< OneD, MultiRegions::ExpListSharedPtr > m_fields
Array holding all dependent variables.
Definition: EquationSystem.h:467
SOLVER_UTILS_EXPORT void SetBoundaryConditions(NekDouble time)
Evaluates the boundary conditions at the given time.
Definition: EquationSystem.cpp:775
std::map< ConstFactorType, NekDouble > ConstFactorMap
Definition: StdRegions.hpp:430
StdRegions::ConstFactorMap factors
Definition: TestVarcoeffHashing.cpp:52
void Smul(int n, const T alpha, const T *x, const int incx, T *y, const int incy)
Scalar multiply y = alpha*x.
Definition: Vmath.hpp:100
References Nektar::StdRegions::eFactorLambda, Nektar::StdRegions::eFactorTau, Nektar::VarcoeffHashingTest::factors, m_epsu, Nektar::SolverUtils::EquationSystem::m_fields, m_varcoeff, Nektar::SolverUtils::EquationSystem::SetBoundaryConditions(), and Vmath::Smul().
Referenced by v_InitObject().
◆ DoOdeRhs()
|
protected |
Computes the reaction terms \(f(u,v)\) and \(g(u,v)\).
Definition at line 231 of file MMFDiffusion.cpp.
234{
236
238 {
240 {
241
242 Array<OneD, NekDouble> x(nq);
243 Array<OneD, NekDouble> y(nq);
245
247
248 for (int k = 0; k < nq; k++)
249 {
253 }
254 }
255 break;
256
258 {
259
260 Array<OneD, NekDouble> x(nq);
261 Array<OneD, NekDouble> y(nq);
263
265
266 for (int k = 0; k < nq; k++)
267 {
268 outarray[0][k] =
272 }
273 }
274 break;
275
277 {
278 Array<OneD, NekDouble> temp(nq);
279
281 NekDouble B = 5.0;
282
283 NekDouble m_a, m_b, m_c, m_d;
284 m_a = B - 1.0;
286 m_c = -1.0 * B;
288
289 temp = Array<OneD, NekDouble>(nq, 0.0);
290 Vmath::Svtvp(nq, m_a, &inarray[0][0], 1, &temp[0], 1, &temp[0], 1);
291 Vmath::Svtvp(nq, m_b, &inarray[1][0], 1, &temp[0], 1,
292 &outarray[0][0], 1);
293
294 temp = Array<OneD, NekDouble>(nq, 0.0);
295 Vmath::Svtvp(nq, m_c, &inarray[0][0], 1, &temp[0], 1, &temp[0], 1);
296 Vmath::Svtvp(nq, m_d, &inarray[1][0], 1, &temp[0], 1,
297 &outarray[1][0], 1);
298 }
299 break;
300
302 {
304 NekDouble B = 5.0;
305
307
308 // cube = phys0*phys0*phy1
309 Array<OneD, NekDouble> cube(nq);
310 Vmath::Vmul(nq, &inarray[0][0], 1, &inarray[0][0], 1, &cube[0], 1);
311 Vmath::Vmul(nq, &inarray[1][0], 1, &cube[0], 1, &cube[0], 1);
312
313 // outarray[0] = A - B*phy0 + phy0*phy0*phy1 - phy0
314 NekDouble coeff = -1.0 * B - 1.0;
315 Array<OneD, NekDouble> tmp(nq);
316 Vmath::Svtvp(nq, coeff, &inarray[0][0], 1, &cube[0], 1, &tmp[0], 1);
317 Vmath::Vadd(nq, &Aonevec[0], 1, &tmp[0], 1, &outarray[0][0], 1);
318
319 // outarray[1] = B*phys0 - phy0*phy0*phy1
320 Vmath::Svtvm(nq, B, &inarray[0][0], 1, &cube[0], 1, &outarray[1][0],
321 1);
322 }
323 break;
324
326 {
327 // \phi - \phi^3/3 - \psi
328 NekDouble a = 0.12;
329 NekDouble b = 0.011;
330 NekDouble c1 = 0.175;
331 NekDouble c2 = 0.03;
333
334 Array<OneD, NekDouble> tmp(nq);
335
336 // Reaction for \phi = c1 \phi ( \phi - a)*(1 - \phi) - c2 v
337 Vmath::Smul(nq, -1.0 * c1, inarray[0], 1, outarray[0], 1);
338 Vmath::Sadd(nq, -1.0 * a, inarray[0], 1, tmp, 1);
339 Vmath::Vmul(nq, tmp, 1, inarray[0], 1, outarray[0], 1);
340 Vmath::Sadd(nq, -1.0, inarray[0], 1, tmp, 1);
341 Vmath::Vmul(nq, tmp, 1, outarray[0], 1, outarray[0], 1);
342
343 Vmath::Smul(nq, -1.0 * c2, inarray[1], 1, tmp, 1);
344 Vmath::Vadd(nq, tmp, 1, outarray[0], 1, outarray[0], 1);
345
346 // Reaction for \psi = b (\phi - d \psi )
348 outarray[1], 1);
349 Vmath::Smul(nq, b, outarray[1], 1, outarray[1], 1);
350 }
351 break;
352
354 {
355 NekDouble a = 0.13;
356 NekDouble b = 0.013;
357 NekDouble c1 = 0.26;
358 NekDouble c2 = 0.1;
360
361 Array<OneD, NekDouble> tmp(nq);
362
363 // Reaction for \phi = c1 \phi ( \phi - a)*(1 - \phi) - c2 u v
364 Vmath::Smul(nq, -1.0 * c1, inarray[0], 1, outarray[0], 1);
365 Vmath::Sadd(nq, -1.0 * a, inarray[0], 1, tmp, 1);
366 Vmath::Vmul(nq, tmp, 1, outarray[0], 1, outarray[0], 1);
367 Vmath::Sadd(nq, -1.0, inarray[0], 1, tmp, 1);
368 Vmath::Vmul(nq, tmp, 1, outarray[0], 1, outarray[0], 1);
369
370 Vmath::Vmul(nq, inarray[0], 1, inarray[1], 1, tmp, 1);
371 Vmath::Smul(nq, -1.0 * c2, tmp, 1, tmp, 1);
372 Vmath::Vadd(nq, tmp, 1, outarray[0], 1, outarray[0], 1);
373
374 // Reaction for \psi = b (\phi - d \psi )
376 outarray[1], 1);
377 Vmath::Smul(nq, b, outarray[1], 1, outarray[1], 1);
378 }
379 break;
380
382 {
383
384 NekDouble a = 0.15;
385 NekDouble c1 = 8.0;
386 NekDouble c2 = 1.0;
387 NekDouble c0 = 0.002;
388 NekDouble mu1 = 0.2;
389 NekDouble mu2 = 0.3;
390
391 Array<OneD, NekDouble> tmp(nq);
392
393 // Reaction for \phi = c1 \phi ( \phi - a)*(1 - \phi) - c2 u v
394 Vmath::Smul(nq, -1.0 * c1, inarray[0], 1, outarray[0], 1);
395 Vmath::Sadd(nq, -1.0 * a, inarray[0], 1, tmp, 1);
396 Vmath::Vmul(nq, tmp, 1, outarray[0], 1, outarray[0], 1);
397 Vmath::Sadd(nq, -1.0, inarray[0], 1, tmp, 1);
398 Vmath::Vmul(nq, tmp, 1, outarray[0], 1, outarray[0], 1);
399
400 Vmath::Vmul(nq, inarray[0], 1, inarray[1], 1, tmp, 1);
401 Vmath::Smul(nq, -1.0 * c2, tmp, 1, tmp, 1);
402 Vmath::Vadd(nq, tmp, 1, outarray[0], 1, outarray[0], 1);
403
404 // Reaction for \psi = (c0 + (\mu1 \psi/(\mu2+\phi) ) )*(-\psi - c1
405 // * \phi*(\phi - a - 1) )
406
407 Vmath::Smul(nq, mu1, inarray[1], 1, outarray[1], 1);
408 Vmath::Sadd(nq, mu2, inarray[0], 1, tmp, 1);
409 Vmath::Vdiv(nq, outarray[1], 1, tmp, 1, outarray[1], 1);
410 Vmath::Sadd(nq, c0, outarray[1], 1, outarray[1], 1);
411
412 Vmath::Sadd(nq, (-a - 1.0), inarray[0], 1, tmp, 1);
413 Vmath::Vmul(nq, inarray[0], 1, tmp, 1, tmp, 1);
414 Vmath::Smul(nq, c1, tmp, 1, tmp, 1);
415 Vmath::Vadd(nq, inarray[1], 1, tmp, 1, tmp, 1);
416 Vmath::Neg(nq, tmp, 1);
417
418 Vmath::Vmul(nq, tmp, 1, outarray[1], 1, outarray[1], 1);
419 }
420 break;
421
422 default:
423 break;
424 }
425}
SOLVER_UTILS_EXPORT int GetTotPoints()
Definition: EquationSystem.h:341
std::vector< double > z(NPUPPER)
std::vector< double > d(NPUPPER *NPUPPER)
void Vmul(int n, const T *x, const int incx, const T *y, const int incy, T *z, const int incz)
Multiply vector z = x*y.
Definition: Vmath.hpp:72
void Svtvp(int n, const T alpha, const T *x, const int incx, const T *y, const int incy, T *z, const int incz)
Svtvp (scalar times vector plus vector): z = alpha*x + y.
Definition: Vmath.hpp:396
void Svtvm(int n, const T alpha, const T *x, const int incx, const T *y, const int incy, T *z, const int incz)
Svtvm (scalar times vector minus vector): z = alpha*x - y.
Definition: Vmath.hpp:424
void Vadd(int n, const T *x, const int incx, const T *y, const int incy, T *z, const int incz)
Add vector z = x+y.
Definition: Vmath.hpp:180
void Vdiv(int n, const T *x, const int incx, const T *y, const int incy, T *z, const int incz)
Multiply vector z = x/y.
Definition: Vmath.hpp:126
void Sadd(int n, const T alpha, const T *x, const int incx, T *y, const int incy)
Add vector y = alpha + x.
Definition: Vmath.hpp:194
References Nektar::UnitTests::d(), Nektar::eFHNAlievPanf, Nektar::eFHNRogers, Nektar::eFHNStandard, Nektar::eTestCube, Nektar::eTestLinearSphere, Nektar::eTestNonlinearSphere, Nektar::eTestPlane, Nektar::SolverUtils::EquationSystem::GetTotPoints(), m_epsilon, Nektar::SolverUtils::EquationSystem::m_fields, Nektar::SolverUtils::MMFSystem::m_pi, m_TestType, Vmath::Neg(), Vmath::Sadd(), Vmath::Smul(), Vmath::Svtvm(), Vmath::Svtvp(), Vmath::Vadd(), Vmath::Vdiv(), Vmath::Vmul(), and Nektar::UnitTests::z().
Referenced by v_InitObject().
◆ Morphogenesis()
|
protected |
Definition at line 542 of file MMFDiffusion.cpp.
544{
546
547 int i, m, n, ind;
548 NekDouble a_n, d_n, gamma_n;
549 NekDouble A_mn, C_mn, theta, phi, radius;
550
551 std::complex<double> Spericharmonic, delta_n, temp;
552 std::complex<double> varphi0, varphi1;
553 std::complex<double> B_mn, D_mn;
554
555 // Set some parameter values
556 int Maxn = 6;
557 int Maxm = 2 * Maxn - 1;
558
560 NekDouble B = 5.0;
561
563 NekDouble m_nu = 0.002;
564
565 NekDouble m_a, m_b, m_c, m_d;
566
567 m_a = B - 1.0;
569 m_c = -1.0 * B;
571
572 Array<OneD, Array<OneD, NekDouble>> Ainit(Maxn);
573 Array<OneD, Array<OneD, NekDouble>> Binit(Maxn);
574
575 for (i = 0; i < Maxn; ++i)
576 {
577 Ainit[i] = Array<OneD, NekDouble>(Maxm, 0.0);
578 Binit[i] = Array<OneD, NekDouble>(Maxm, 0.0);
579 }
580
581 Ainit[5][0] = -0.5839;
582 Ainit[5][1] = -0.8436;
583 Ainit[5][2] = -0.4764;
584 Ainit[5][3] = 0.6475;
585 Ainit[5][4] = 0.1886;
586 Ainit[5][5] = 0.8709;
587 Ainit[5][6] = -0.8338;
588 Ainit[5][7] = 0.1795;
589 Ainit[5][8] = -0.7873;
590 Ainit[5][9] = 0.8842;
591 Ainit[5][10] = 0.2943;
592
593 Binit[5][0] = -0.6263;
594 Binit[5][1] = 0.9803;
595 Binit[5][2] = 0.7222;
596 Binit[5][3] = 0.5945;
597 Binit[5][4] = 0.6026;
598 Binit[5][5] = -0.2076;
599 Binit[5][6] = 0.4556;
600 Binit[5][7] = 0.6024;
601 Binit[5][8] = 0.9695;
602 Binit[5][9] = -0.4936;
603 Binit[5][10] = 0.1098;
604
605 Array<OneD, NekDouble> u(nq);
606 Array<OneD, NekDouble> v(nq);
607 Array<OneD, NekDouble> x(nq);
608 Array<OneD, NekDouble> y(nq);
610
612 for (int i = 0; i < nq; ++i)
613 {
615
616 // theta is in [0, pi]
618
619 // phi is in [0, 2*pi]
620 phi = atan2(y[i], x[i]) + m_pi;
621
622 varphi0 = 0.0 * varphi0;
623 varphi1 = 0.0 * varphi1;
624 for (n = 0; n < Maxn; ++n)
625 {
626 // Set up parameters
627 a_n = m_a - m_mu * (n * (n + 1) / radius / radius);
628 d_n = m_d - m_nu * (n * (n + 1) / radius / radius);
629
630 gamma_n = 0.5 * (a_n + d_n);
631
632 temp = (a_n + d_n) * (a_n + d_n) - 4.0 * (a_n * d_n - m_b * m_c);
633 delta_n = 0.5 * sqrt(temp);
634
635 for (m = -n; m <= n; ++m)
636 {
637 ind = m + n;
638 A_mn = Ainit[n][ind];
639 C_mn = Binit[n][ind];
640
641 B_mn = ((a_n - gamma_n) * Ainit[n][ind] + m_b * Binit[n][ind]) /
642 delta_n;
643 D_mn = (m_c * Ainit[n][ind] + (d_n - gamma_n) * Binit[n][ind]) /
644 delta_n;
645
646 Spericharmonic =
647 boost::math::spherical_harmonic(n, m, theta, phi);
648 varphi0 += exp(gamma_n * time) *
649 (A_mn * cosh(delta_n * time) +
650 B_mn * sinh(delta_n * time)) *
651 Spericharmonic;
652 varphi1 += exp(gamma_n * time) *
653 (C_mn * cosh(delta_n * time) +
654 D_mn * sinh(delta_n * time)) *
655 Spericharmonic;
656 }
657 }
658
659 u[i] = varphi0.real();
660 v[i] = varphi1.real();
661 }
662
664 {
665 case 0:
666 {
667 outfield = u;
668 }
669 break;
670
671 case 1:
672 {
673 outfield = v;
674 }
675 break;
676 }
677}
References FilterPython_Function::field, Nektar::SolverUtils::EquationSystem::GetTotPoints(), Nektar::SolverUtils::EquationSystem::m_fields, m_mu, Nektar::SolverUtils::MMFSystem::m_pi, tinysimd::sqrt(), and Nektar::UnitTests::z().
Referenced by v_EvaluateExactSolution(), and v_SetInitialConditions().
◆ PlanePhiWave()
Definition at line 679 of file MMFDiffusion.cpp.
680{
682 Array<OneD, NekDouble> outarray(nq, 0.0);
683
684 Array<OneD, NekDouble> x(nq);
685 Array<OneD, NekDouble> y(nq);
687
689
690 NekDouble xmin, ymin, xmax;
691
692 xmin = Vmath::Vmin(nq, x, 1);
693 xmax = Vmath::Vmax(nq, x, 1);
694 ymin = Vmath::Vmin(nq, y, 1);
695
696 NekDouble xp, yp, xp2;
697 for (int i = 0; i < nq; i++)
698 {
700 {
702 {
703 NekDouble radiusofinit = 4.0;
704 NekDouble frontstiff = 0.1;
705
706 xp = x[i] - xmin;
707 outarray[i] =
708 1.0 / (1.0 + exp((xp - radiusofinit) / frontstiff));
709 }
710 break;
711
713 {
714 NekDouble radiusofinit = 3.0;
715 NekDouble frontstiff = 0.1;
716
717 xp = x[i] - xmin;
718 xp2 = x[i] - xmax;
719
720 outarray[i] =
721 1.0 / (1.0 +
722 exp((sqrt(xp * xp) - radiusofinit) / frontstiff)) +
723 1.0 / (1.0 +
724 exp((sqrt(xp2 * xp2) - radiusofinit) / frontstiff));
725 }
726 break;
727
729 {
730 NekDouble radiusofinit = 6.0;
731 NekDouble frontstiff = 0.1;
732
733 xp = x[i] - xmin;
734 outarray[i] =
735 1.0 / (1.0 + exp((xp - radiusofinit) / frontstiff));
736 }
737 break;
738
740 {
741 NekDouble radiusofinit = 6.0;
742 NekDouble frontstiff = 0.1;
743 NekDouble bs = 2.0;
744
745 xp = x[i] - xmin;
746 yp = y[i] - ymin;
747 outarray[i] =
748 1.0 /
749 (1.0 + exp((sqrt(xp * xp + yp * yp) / bs - radiusofinit) /
750 frontstiff));
751 }
752 break;
753
755 {
757 NekDouble radiusofinit = 10.0;
758
759 xloc = x[i] - m_InitPtx;
760 yloc = y[i] - m_InitPty;
762
764
765 xloc = xloc / radiusofinit;
766 yloc = yloc / radiusofinit;
767 zloc = zloc / radiusofinit;
768
770 {
771 outarray[i] =
772 exp(-(1.0 / 2.0) *
773 (xloc * xloc + yloc * yloc + zloc * zloc));
774 }
775
776 else
777 {
778 outarray[i] = 0.0;
779 }
780 }
781 break;
782
784 {
785 NekDouble radiusofinit = 3.0;
786 NekDouble frontstiff = 0.1;
787 xp = x[i] - 4.0;
788 yp = y[i];
789 outarray[i] =
790 (1.0 / (1.0 + exp(2.0 * yp))) *
791 (1.0 / (1.0 + exp(-2.0 * xp))) *
792 (1.0 / (1.0 + exp((xp - radiusofinit) / frontstiff)));
793 }
794 break;
795
796 default:
797 break;
798 }
799 }
800
801 return outarray;
802}
static NekDouble rad(NekDouble x, NekDouble y)
Definition: Interpreter/Interpreter.cpp:76
T Vmin(int n, const T *x, const int incx)
Return the minimum element in x - called vmin to avoid conflict with min.
Definition: Vmath.hpp:725
T Vmax(int n, const T *x, const int incx)
Return the maximum element in x – called vmax to avoid conflict with max.
Definition: Vmath.hpp:644
References Nektar::eBothEnds, Nektar::eCenter, Nektar::eLeft, Nektar::eLeftBottomCorner, Nektar::ePoint, Nektar::eSpiralDock, Nektar::SolverUtils::EquationSystem::GetTotPoints(), Nektar::SolverUtils::EquationSystem::m_fields, m_InitPtx, m_InitPty, m_InitPtz, m_InitWaveType, Nektar::LibUtilities::rad(), tinysimd::sqrt(), Vmath::Vmax(), Vmath::Vmin(), and Nektar::UnitTests::z().
Referenced by v_SetInitialConditions().
◆ TestCubeProblem()
|
protected |
Definition at line 522 of file MMFDiffusion.cpp.
525{
527
528 Array<OneD, NekDouble> x(nq);
529 Array<OneD, NekDouble> y(nq);
531
533
534 outfield = Array<OneD, NekDouble>(nq);
535 for (int k = 0; k < nq; k++)
536 {
539 }
540}
References Nektar::SolverUtils::EquationSystem::GetTotPoints(), Nektar::SolverUtils::EquationSystem::m_fields, Nektar::SolverUtils::MMFSystem::m_pi, and Nektar::UnitTests::z().
Referenced by v_EvaluateExactSolution(), and v_SetInitialConditions().
◆ TestPlaneProblem()
|
protected |
Definition at line 502 of file MMFDiffusion.cpp.
505{
507
508 Array<OneD, NekDouble> x(nq);
509 Array<OneD, NekDouble> y(nq);
511
513
514 outfield = Array<OneD, NekDouble>(nq);
515 for (int k = 0; k < nq; k++)
516 {
518 cos(m_pi * y[k]);
519 }
520}
References Nektar::SolverUtils::EquationSystem::GetTotPoints(), Nektar::SolverUtils::EquationSystem::m_fields, Nektar::SolverUtils::MMFSystem::m_pi, and Nektar::UnitTests::z().
Referenced by v_EvaluateExactSolution(), and v_SetInitialConditions().
◆ v_EvaluateExactSolution()
|
overrideprotectedvirtual |
Reimplemented from Nektar::SolverUtils::EquationSystem.
Definition at line 804 of file MMFDiffusion.cpp.
807{
809 {
811 {
812 TestPlaneProblem(time, outfield);
813 }
814 break;
815
817 {
818 TestCubeProblem(time, outfield);
819 }
820 break;
821
824 {
826 }
827 break;
828
832 {
834 outfield = Array<OneD, NekDouble>(nq, 0.0);
835 }
836 /* Falls through. */
837 default:
838 {
840 }
841 break;
842 }
843}
void Morphogenesis(const NekDouble time, unsigned int field, Array< OneD, NekDouble > &outfield)
Definition: MMFDiffusion.cpp:542
void TestCubeProblem(const NekDouble time, Array< OneD, NekDouble > &outfield)
Definition: MMFDiffusion.cpp:522
void TestPlaneProblem(const NekDouble time, Array< OneD, NekDouble > &outfield)
Definition: MMFDiffusion.cpp:502
virtual SOLVER_UTILS_EXPORT void v_EvaluateExactSolution(unsigned int field, Array< OneD, NekDouble > &outfield, const NekDouble time)
Definition: EquationSystem.cpp:1141
References Nektar::eFHNAlievPanf, Nektar::eFHNRogers, Nektar::eFHNStandard, Nektar::eTestCube, Nektar::eTestLinearSphere, Nektar::eTestNonlinearSphere, Nektar::eTestPlane, FilterPython_Function::field, Nektar::SolverUtils::EquationSystem::GetTotPoints(), m_TestType, Morphogenesis(), TestCubeProblem(), TestPlaneProblem(), and Nektar::SolverUtils::EquationSystem::v_EvaluateExactSolution().
◆ v_GenerateSummary()
|
overrideprotectedvirtual |
Prints a summary of the model parameters.
Reimplemented from Nektar::SolverUtils::MMFSystem.
Definition at line 845 of file MMFDiffusion.cpp.
846{
847 MMFSystem::v_GenerateSummary(s);
853 {
856 }
857}
SOLVER_UTILS_EXPORT void v_GenerateSummary(SummaryList &s) override
Virtual function for generating summary information.
Definition: MMFSystem.cpp:2463
void AddSummaryItem(SummaryList &l, const std::string &name, const std::string &value)
Adds a summary item to the summary info list.
Definition: Misc.cpp:47
References Nektar::SolverUtils::AddSummaryItem(), Nektar::eTestLinearSphere, m_epsilon, m_epsu, m_TestType, Nektar::TestTypeMap, and Nektar::SolverUtils::MMFSystem::v_GenerateSummary().
◆ v_InitObject()
|
overrideprotectedvirtual |
Init object for UnsteadySystem class.
Initialization object for UnsteadySystem class.
Reimplemented from Nektar::SolverUtils::UnsteadySystem.
Definition at line 64 of file MMFDiffusion.cpp.
65{
66 UnsteadySystem::v_InitObject(DeclareFields);
67
70 int MFdim = 3;
71
72 // Diffusivity coefficient for e^j
77
78 // Diffusivity coefficient for u^j
82
86
88 Array<OneD, Array<OneD, NekDouble>> Anisotropy(shapedim);
89 for (int j = 0; j < shapedim; ++j)
90 {
91 Anisotropy[j] = Array<OneD, NekDouble>(nq, 1.0);
93 }
94
95 MMFSystem::MMFInitObject(Anisotropy);
96
97 // Define ProblemType
99 {
101 int i;
103 {
105 {
107 break;
108 }
109 }
110 }
111 else
112 {
114 }
115
117 {
120 {
122 {
124 break;
125 }
126 }
127 }
128 else
129 {
131 }
132
133 StdRegions::VarCoeffType MMFCoeffs[15] = {
142
143 int indx;
144 Array<OneD, NekDouble> tmp(nq);
145 for (int k = 0; k < MFdim; ++k)
146 {
147 // For Moving Frames
148 indx = 5 * k;
149
151 {
153 m_varcoeff[MMFCoeffs[indx + j]] = tmp;
154 }
155
156 // m_DivMF
158 m_varcoeff[MMFCoeffs[indx + 3]] = tmp;
159
160 // \| e^k \|
161 tmp = Array<OneD, NekDouble>(nq, 0.0);
163 {
165 &m_movingframes[k][i * nq], 1, &tmp[0], 1, &tmp[0], 1);
166 }
167
168 m_varcoeff[MMFCoeffs[indx + 4]] = tmp;
169 }
170
172 {
174 }
175
177}
void DefineOdeRhs(FuncPointerT func, ObjectPointerT obj)
Definition: TimeIntegrationSchemeOperators.h:84
void DefineImplicitSolve(FuncPointerT func, ObjectPointerT obj)
Definition: TimeIntegrationSchemeOperators.h:100
void DoOdeRhs(const Array< OneD, const Array< OneD, NekDouble > > &inarray, Array< OneD, Array< OneD, NekDouble > > &outarray, const NekDouble time)
Computes the reaction terms and .
Definition: MMFDiffusion.cpp:231
void DoImplicitSolve(const Array< OneD, const Array< OneD, NekDouble > > &inarray, Array< OneD, Array< OneD, NekDouble > > &outarray, NekDouble time, NekDouble lambda)
Solve for the diffusion term.
Definition: MMFDiffusion.cpp:185
Array< OneD, Array< OneD, NekDouble > > m_DivMF
Definition: MMFSystem.h:192
Array< OneD, Array< OneD, NekDouble > > m_movingframes
Definition: MMFSystem.h:183
SOLVER_UTILS_EXPORT void MMFInitObject(const Array< OneD, const Array< OneD, NekDouble > > &Anisotropy, const int TangentXelem=-1)
Definition: MMFSystem.cpp:51
LibUtilities::TimeIntegrationSchemeOperators m_ode
The time integration scheme operators to use.
Definition: UnsteadySystem.h:89
bool m_explicitDiffusion
Indicates if explicit or implicit treatment of diffusion is used.
Definition: UnsteadySystem.h:107
SOLVER_UTILS_EXPORT void v_InitObject(bool DeclareField=true) override
Init object for UnsteadySystem class.
Definition: UnsteadySystem.cpp:85
void Vvtvp(int n, const T *w, const int incw, const T *x, const int incx, const T *y, const int incy, T *z, const int incz)
vvtvp (vector times vector plus vector): z = w*x + y
Definition: Vmath.hpp:366
void Fill(int n, const T alpha, T *x, const int incx)
Fill a vector with a constant value.
Definition: Vmath.hpp:54
void Vcopy(int n, const T *x, const int incx, T *y, const int incy)
Definition: Vmath.hpp:825
References Nektar::LibUtilities::TimeIntegrationSchemeOperators::DefineImplicitSolve(), Nektar::LibUtilities::TimeIntegrationSchemeOperators::DefineOdeRhs(), DoImplicitSolve(), DoOdeRhs(), Nektar::StdRegions::eVarCoeffMF1Div, Nektar::StdRegions::eVarCoeffMF1Mag, Nektar::StdRegions::eVarCoeffMF1x, Nektar::StdRegions::eVarCoeffMF1y, Nektar::StdRegions::eVarCoeffMF1z, Nektar::StdRegions::eVarCoeffMF2Div, Nektar::StdRegions::eVarCoeffMF2Mag, Nektar::StdRegions::eVarCoeffMF2x, Nektar::StdRegions::eVarCoeffMF2y, Nektar::StdRegions::eVarCoeffMF2z, Nektar::StdRegions::eVarCoeffMF3Div, Nektar::StdRegions::eVarCoeffMF3Mag, Nektar::StdRegions::eVarCoeffMF3x, Nektar::StdRegions::eVarCoeffMF3y, Nektar::StdRegions::eVarCoeffMF3z, Vmath::Fill(), Nektar::InitWaveTypeMap, Nektar::SolverUtils::MMFSystem::m_DivMF, m_epsilon, m_epsu, Nektar::SolverUtils::UnsteadySystem::m_explicitDiffusion, Nektar::SolverUtils::EquationSystem::m_fields, m_InitPtx, m_InitPty, m_InitPtz, m_InitWaveType, Nektar::SolverUtils::MMFSystem::m_movingframes, Nektar::SolverUtils::UnsteadySystem::m_ode, Nektar::SolverUtils::EquationSystem::m_session, Nektar::SolverUtils::EquationSystem::m_spacedim, m_TestType, m_varcoeff, Nektar::SolverUtils::MMFSystem::MMFInitObject(), Nektar::SIZE_TestType, tinysimd::sqrt(), Nektar::TestTypeMap, Nektar::SolverUtils::UnsteadySystem::v_InitObject(), Vmath::Vcopy(), and Vmath::Vvtvp().
◆ v_SetInitialConditions()
|
overrideprotectedvirtual |
Sets a custom initial condition.
Reimplemented from Nektar::SolverUtils::EquationSystem.
Definition at line 430 of file MMFDiffusion.cpp.
433{
435
437 {
439 {
440 Array<OneD, NekDouble> u(nq);
441
442 TestPlaneProblem(initialtime, u);
443 m_fields[0]->SetPhys(u);
444 }
445 break;
446
448 {
449 Array<OneD, NekDouble> u(nq);
450
451 TestCubeProblem(initialtime, u);
452 m_fields[0]->SetPhys(u);
453 }
454 break;
455
458 {
459 Array<OneD, NekDouble> u(nq);
460 Array<OneD, NekDouble> v(nq);
461
462 Morphogenesis(initialtime, 0, u);
463 Morphogenesis(initialtime, 1, v);
464
465 m_fields[0]->SetPhys(u);
466 m_fields[1]->SetPhys(v);
467 }
468 break;
469
473 {
477 }
478 break;
479
480 default:
481 {
483 }
484 break;
485 }
486
487 // forward transform to fill the modal coeffs
489 {
492 m_fields[i]->UpdateCoeffs());
493 }
494
495 if (dumpInitialConditions)
496 {
498 WriteFld(outname);
499 }
500}
Array< OneD, NekDouble > PlanePhiWave()
Definition: MMFDiffusion.cpp:679
virtual SOLVER_UTILS_EXPORT void v_SetInitialConditions(NekDouble initialtime=0.0, bool dumpInitialConditions=true, const int domain=0)
Definition: EquationSystem.cpp:1063
SOLVER_UTILS_EXPORT void WriteFld(const std::string &outname)
Write field data to the given filename.
Definition: EquationSystem.cpp:1315
References Nektar::eFHNAlievPanf, Nektar::eFHNRogers, Nektar::eFHNStandard, Nektar::eTestCube, Nektar::eTestLinearSphere, Nektar::eTestNonlinearSphere, Nektar::eTestPlane, Nektar::SolverUtils::EquationSystem::GetTotPoints(), Nektar::SolverUtils::EquationSystem::m_fields, Nektar::SolverUtils::EquationSystem::m_sessionName, m_TestType, Morphogenesis(), PlanePhiWave(), TestCubeProblem(), TestPlaneProblem(), Nektar::SolverUtils::EquationSystem::v_SetInitialConditions(), Nektar::SolverUtils::EquationSystem::WriteFld(), and Vmath::Zero().
Friends And Related Function Documentation
◆ MemoryManager< MMFDiffusion >
|
friend |
Definition at line 74 of file MMFDiffusion.h.
Member Data Documentation
◆ className
|
static |
Initial value:
=
tKey RegisterCreatorFunction(tKey idKey, CreatorFunction classCreator, std::string pDesc="")
Register a class with the factory.
Definition: BasicUtils/NekFactory.hpp:197
static SolverUtils::EquationSystemSharedPtr create(const LibUtilities::SessionReaderSharedPtr &pSession, const SpatialDomains::MeshGraphSharedPtr &pGraph)
Creates an instance of this class.
Definition: MMFDiffusion.h:85
EquationSystemFactory & GetEquationSystemFactory()
Definition: EquationSystem.cpp:99
Name of class.
Definition at line 96 of file MMFDiffusion.h.
◆ m_epsilon
Definition at line 148 of file MMFDiffusion.h.
Referenced by DoOdeRhs(), v_GenerateSummary(), and v_InitObject().
◆ m_epsu
Definition at line 149 of file MMFDiffusion.h.
Referenced by DoImplicitSolve(), v_GenerateSummary(), and v_InitObject().
◆ m_InitPtx
|
protected |
Definition at line 101 of file MMFDiffusion.h.
Referenced by PlanePhiWave(), and v_InitObject().
◆ m_InitPty
|
protected |
Definition at line 101 of file MMFDiffusion.h.
Referenced by PlanePhiWave(), and v_InitObject().
◆ m_InitPtz
|
protected |
Definition at line 101 of file MMFDiffusion.h.
Referenced by PlanePhiWave(), and v_InitObject().
◆ m_InitWaveType
|
protected |
Definition at line 100 of file MMFDiffusion.h.
Referenced by PlanePhiWave(), and v_InitObject().
◆ m_TestType
|
protected |
Definition at line 99 of file MMFDiffusion.h.
Referenced by DoOdeRhs(), v_EvaluateExactSolution(), v_GenerateSummary(), v_InitObject(), and v_SetInitialConditions().
◆ m_varcoeff
|
private |
Variable diffusivity.
Definition at line 147 of file MMFDiffusion.h.
Referenced by DoImplicitSolve(), and v_InitObject().
