PviewScripts/pyAlya_MpioAvgTurb.py at master · kumarvrai/PviewScripts · GitHub

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
#!/bin/env python
#
# HiFiTurb database computation.
#
# Last rev: 19/02/2021
from __future__ import print_function, division

# Please do not delete this part otherwise it will not work
# you have been warned after a long weekend of debugging
import mpi4py
mpi4py.rc.recv_mprobe = False

import numpy as np
import pyAlya
import os
import sys


# Parameters
rho, mu = 1.0, float(sys.argv[5])
lam     = 0.01

BASEDIR        = ''
ALT_BASEDIR    = ''
CASESTR        = sys.argv[1]
VARLIST        = ['AVVEL','AVPRE']
START, DT, END = int(sys.argv[2]),int(sys.argv[3]),int(sys.argv[4])

FILE_FMT = sys.argv[6]
COMM = sys.argv[7]
PLOTVAR = sys.argv[8]

#Create list of instances
if(START==END):
  listOfInstants = [END]
else:
  listOfInstants = [ii for ii in range(START,END+DT,DT)]


## Create the subdomain mesh
mesh = pyAlya.Mesh.read(CASESTR,basedir=BASEDIR,alt_basedir=ALT_BASEDIR,fmt=FILE_FMT,read_commu=True if COMM == 1 else False,read_massm=False)


pyAlya.pprint(0,'Run (%d instants)...' % len(listOfInstants),flush=True)

time=0
for instant in listOfInstants:
  pyAlya.pprint(1,'First loop, instant %d...'%instant,flush=True)
  # Read field
  field, header = pyAlya.Field.read(CASESTR,VARLIST,instant,mesh.xyz,basedir=BASEDIR,fmt=FILE_FMT)

  # Compute and smooth the gradient of velocity
  #gradv = mesh.smooth(mesh.gradient(field['AVVEL']),iters=3)
  gradv = mesh.gradient(field['AVVEL'])

  # # Store gradients
  # field['GRADV'] = mesh.newArray(ndim=6)
  # field['GRADV'][:,0] = gradv[:,0] # XX
  # field['GRADV'][:,1] = gradv[:,4] # YY
  # field['GRADV'][:,2] = gradv[:,8] # ZZ
  # field['GRADV'][:,3] = gradv[:,1] # XY
  # field['GRADV'][:,4] = gradv[:,2] # XZ
  # field['GRADV'][:,5] = gradv[:,5] # YZ

  # Compute Vorticity, Q and Omega from the gradient
  if('ALL' not in PLOTVAR):
    if('VOR' in PLOTVAR):
      field['VORTI'] = pyAlya.postproc.vorticity(gradv)
    if('QCR' in PLOTVAR):
      field['QCRIT'] = pyAlya.postproc.QCriterion(gradv)
    if('LAM' in PLOTVAR):
      field['LAMB2'] = pyAlya.postproc.Lambda2Criterion(gradv)
    if('OMG' in PLOTVAR):
      field['OMEGA'] = pyAlya.postproc.OmegaCriterion(gradv,epsilon=0.001,modified=False)
    if('ROR' in PLOTVAR):
      field['RORTX'] = pyAlya.postproc.RortexCriterion(gradv)
    if('ORX' in PLOTVAR):
      field['OMERX'] = pyAlya.postproc.OmegaRortexCriterion(gradv,epsilon=0.001,modified=False)
  else:
    field['VORTI'] = pyAlya.postproc.vorticity(gradv)
    field['QCRIT'] = pyAlya.postproc.QCriterion(gradv)
    field['LAMB2'] = pyAlya.postproc.Lambda2Criterion(gradv)
    field['OMEGA'] = pyAlya.postproc.OmegaCriterion(gradv,epsilon=0.001,modified=False)
    field['RORTX'] = pyAlya.postproc.RortexCriterion(gradv)
    field['OMERX'] = pyAlya.postproc.OmegaRortexCriterion(gradv,epsilon=0.001,modified=False)

  #Write the fields
  pyAlya.pprint(1,'Writing MPIO...',flush=True)
  field.write(CASESTR,instant,header.time,basedir=BASEDIR,fmt=FILE_FMT,exclude_vars=VARLIST)

pyAlya.cr_info()