MultirateTsunami.py 3.29 KB
Newer Older
bseny's avatar
bseny committed
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
from dgpy import *
from PartitionTsunami import *
from InputTsunami import *
from math import *
import time, os, sys



model = GModel()
loadAndPartitionMesh(model, filename, mL, algo)

if(not os.path.exists(bathname)):
  Msg.Fatal('No valid bathymetry file :' + bathname + ' found, please launch "rundgpy DiffuseBathymetryTsunami.py"')

stereo = dgSpaceTransformSpherical(6371220.)
groups = dgGroupCollection(model, dimension, order, stereo)
XYZ = groups.getFunctionCoordinates();
claw = dgConservationLawShallowWater2d()
bath = dgDofContainer(groups, 1);
20
bath.importIdx(bathname)
bseny's avatar
bseny committed
21
solution = dgDofContainer(groups, claw.getNbFields())
bseny's avatar
bseny committed
22
Finit = functionC(tlib, "initial_condition_okada_honshu",3, [XYZ])
bseny's avatar
bseny committed
23
solution.interpolate(Finit)
bseny's avatar
bseny committed
24
FCoriolis= functionC(tlib, "coriolis", 1, [XYZ])
bseny's avatar
bseny committed
25
26
27
28
29
30
31
32
33
claw.setBathymetry(bath.getFunction())
claw.setBathymetryGradient(bath.getFunctionGradient())
claw.setCoriolisFactor(FCoriolis)
claw.setIsLinear(True)
claw.setIsSpherical(1, XYZ)
wall=claw.newBoundaryWall()
claw.addBoundaryCondition("Coast", wall)


bseny's avatar
bseny committed
34
bottomDrag = functionC(tlib, "bottomDrag", 1, [function.getSolution(), bath.getFunction()])
bseny's avatar
bseny committed
35
36
claw.setLinearDissipation(bottomDrag)

bseny's avatar
bseny committed
37
38
39
rk = dgMultirateERK(groups, claw, RK_TYPE)
dt=rk.splitGroupsForMultirate(mL, solution, [solution, bath])
rk.printMultirateInfo(0)
bseny's avatar
bseny committed
40
41
42
43
44
45

evalSol=dgFunctionEvaluator(groups, solution.getFunction())

nbExport = 0
nbSteps=int(ceil((Tf-Ti)/dt))

bseny's avatar
bseny committed
46
exporterSol=dgIdxExporter(solution, outputDir+'/mr-tsunami')
bseny's avatar
bseny committed
47
48
49
50
51
52
53
norm=solution.norm()

if (Msg.GetCommRank() == 0):
  print ''
  print(colored('******************** Initializing ********************', "blue"))
  print ''
  print '     - Number of processors:', Msg.GetCommSize()
bseny's avatar
bseny committed
54
  print '     - Multirate Runge Kutta Scheme:', RK_TYPE
bseny's avatar
bseny committed
55
56
57
58
  print '     - Initial Time (Ti):',  printTime(Ti)
  print '     - Final Time (Tf):',  printTime(Tf)
  print '     - Reference Time Step:', printTime(dt)
  print '     - Number of Iterartions:', nbSteps
59
  print '     - Theoretical Speedup=', rk.speedUp()
bseny's avatar
bseny committed
60
61
62
63
64
65
66
67
  print '     - Norm of Solution at Ti:',  norm
  print ''
  print(colored('******************************************************', "blue"))
  print ''

t=Ti
exporterSol.exportIdx(0, t)
startcpu = time.time()
bseny's avatar
bseny committed
68
69
70
71
for i in range (0,nbSteps) :
  if(i == nbSteps - 1):
    dt = Tf -t
  rk.iterate(solution, dt, t)
bseny's avatar
bseny committed
72
73
74
75
76
77
78
79
80
81
82
83
84
  t = t + dt
  if ( i%iter == 0 ):
     norm=solution.norm()
     if(Msg.GetCommRank() == 0):
       print(colored('|ITER| %d',  "red")%(i))
       print(colored('------------------------------------------------------', "red"))
       print'|TIME|', printTime(t)
       print'|NORM|',norm
       print'|CPUT|',printTime(time.time() - startcpu)
       print(colored('------------------------------------------------------', "red"))
       print''
  if ( i%export  == 0 ):
    exporterSol.exportIdx(i, t)
bseny's avatar
bseny committed
85
    #solution.exportFunctionMsh(currentFunction, outputDir+'/mr_tsunami_current-%06d'%(i), t, nbExport, '')
bseny's avatar
bseny committed
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
    nbExport = nbExport + 1

Msg.Barrier()
endcpu=time.time()

if (Msg.GetCommRank() == 0):
  print ''
  print(colored('********************     End      ********************', "blue"))
  print ''
  print '     - Final Time (Tf):',  printTime(t)
  print '     - Norm of Solution at Tf:',  norm
  print '     - Total CPU Time:', endcpu-startcpu
  print ''
  print(colored('******************************************************', "blue"))
  print ''
Msg.Exit(0)