No commit message

No commit message

HySoP/hysop/operator/tests/test_diff_poisson.py

-# -*- coding: utf-8 -*-
-import time
-
-import parmepy as pp
-import numpy as np
-from parmepy.constants import PARMES_REAL, ORDER
-from parmepy.fields.analytical import AnalyticalField
-from parmepy.operator.poisson import Poisson
-from parmepy.operator.diffusion import Diffusion
-from parmepy.problem.navier_stokes import NSProblem
-from parmepy.operator.monitors.energy_enstrophy import Energy_enstrophy
-from math import sqrt, pi, exp
-
-
-def computeVel(x, y, z):
-    vx = 0.
-    vy = 0.
-    vz = 0.
-    return vx, vy, vz
-
-def computeVort(x, y, z):
-    xc = 1. / 2.
-    yc = 1. / 2.
-    zc = 1. / 4.
-    R = 0.2
-    sigma = R / 2.
-    Gamma = 0.0075
-    dist = sqrt((x-xc) ** 2 + (y-yc) ** 2)
-    s2 = (z - zc) ** 2 + (dist - R) ** 2
-    wx = 0.
-    wy = 0.
-    wz = 0.
-    if (dist != 0.):
-        cosTheta = (x - xc) / dist
-        sinTheta = (y - yc) / dist
-        wTheta = Gamma / (pi * sigma ** 2) * \
-                 exp(-(s2 / sigma ** 2))
-        wx = - wTheta * sinTheta
-        wy = wTheta * cosTheta
-        wz = 0.
-    return wx, wy, wz
-
-def test_Diff_Poisson():
-    # Parameters
-    nb = 65
-    dim = 3
-    boxLength = [1., 1., 1.]
-    boxMin = [0., 0., 0.]
-    nbElem = [nb, nb, nb]
-
-    timeStep = 0.01
-    finalTime = 150 * timeStep
-    outputFilePrefix = './Energies'
-    outputModulo = 10
-
-    t0 = time.time()
-
-    ## Domain
-    box = pp.Box(dim, length=boxLength, origin=boxMin)
-
-    ## Fields
-    velo = AnalyticalField(domain=box, formula=computeVel, 
-                           name='Velocity', vector=True)
-    vorti = AnalyticalField(domain=box, formula=computeVort, 
-                            name='Vorticity', vector=True)
-
-    ## FFT Diffusion operators and FFT Poisson solver
-    diffusion = Diffusion(velo, vorti,
-                          resolutions={velo: nbElem,
-                                       vorti: nbElem},
-                          method='',
-                          viscosity=0.002, 
-                          with_curl=False
-                         )
-
-    poisson = Poisson(velo, vorti,
-                      resolutions={velo: nbElem,
-                                   vorti: nbElem},
-                      method='',
-                      domain=box
-                     )
-
-    ## Problem
-    pb = NSProblem([diffusion, poisson],
-                    monitors=[Energy_enstrophy(
-                                fields=[velo, vorti],
-                                frequency=outputModulo,
-                                outputPrefix=outputFilePrefix)])
-
-
-    ## Setting solver to Problem
-    pb.setUp(finalTime, timeStep)
-
-    t1 = time.time()
-
-    ## Solve problem
-    timings = pb.solve()
-
-    tf = time.time()
-
-    print "\n"
-    print "Total time : ", tf - t0, "sec (CPU)"
-    print "Init time : ", t1 - t0, "sec (CPU)"
-    print "Solving time : ", tf - t1, "sec (CPU)"
-
-
-if __name__ == "__main__":
-    test_Diff_Poisson()