erdc / proteus

A computational methods and simulation toolkit
http://proteustoolkit.org
MIT License
88 stars 56 forks source link

TPF: initial conditions refactoring (ctd.) #1194

Closed zhang-alvin closed 4 years ago

zhang-alvin commented 4 years ago

Mandatory Checklist

Please ensure that the following criteria are met:

As a general rule of thumb, try to follow PEP8 guidelines.

Description

See #1071 for initial PR.

Sample Dambreak case using RAN2P, VOF, NCLS, RD, MCorr

"""
dambreak 2-D
"""
from proteus import (Domain, Context)
from proteus.mprans.SpatialTools import Tank2D
from proteus.mprans import SpatialTools as st
import proteus.TwoPhaseFlow.TwoPhaseFlowProblem as TpFlow
from proteus.Gauges import PointGauges, LineIntegralGauges, LineGauges
import numpy as np

# *************************** #
# ***** GENERAL OPTIONS ***** #
# *************************** #
opts= Context.Options([
    ("final_time",2.0,"Final time for simulation"),
    ("dt_output",0.01,"Time interval to output solution"),
    ("cfl",0.9,"Desired CFL restriction"),
    ("he",0.08,"he relative to Length of domain in x"),
    ("x_tank",3.22,"extent of domain in x"),
    ("y_tank",1.8,"extent of domain in y")
    ])

#Try adding these variables as Context options#
waterLine_y = 0.6 #the extent of the water column in +y from 0
waterLine_x = 1.2 #the extent of the water column in +x from 0

# *************************** #
# ***** DOMAIN AND MESH ***** #
# ****************** #******* #
tank_dim = (opts.x_tank,opts.y_tank)

structured=False
if structured:
    domain = Domain.RectangularDomain(tank_dim)
else:
    domain = Domain.PlanarStraightLineGraphDomain()

# ----- TANK ----- #
tank = Tank2D(domain, tank_dim)

# ----- BOUNDARY CONDITIONS ----- #
tank.BC['y+'].setAtmosphere()
tank.BC['y-'].setFreeSlip()
tank.BC['x+'].setFreeSlip()
tank.BC['x-'].setFreeSlip()

# ****************************** #
# ***** INITIAL CONDITIONS ***** #
# ****************************** #
class zero(object):
    def uOfXT(self,x,t):
        return 0.0

waterLine_y = 0.6
waterLine_x = 1.2

class PHI_IC:
    def uOfXT(self, x, t):
        phi_x = x[0] - waterLine_x
        phi_y = x[1] - waterLine_y
        if phi_x < 0.0:
            if phi_y < 0.0:
                return max(phi_x, phi_y)
            else:
                return phi_y
        else:
            if phi_y < 0.0:
                return phi_x
            else:
                return (phi_x ** 2 + phi_y ** 2)**0.5

class VF_IC:
    def __init__(self):
        self.phi = PHI_IC()
    def uOfXT(self, x, t):
        from proteus.ctransportCoefficients import smoothedHeaviside
        return smoothedHeaviside(1.5*opts.he,self.phi.uOfXT(x,0.0))

#########################
# ***** Numerics ****** #
#########################
#for structured
nny = int(tank_dim[1]/opts.he)+1
nnx = int(tank_dim[0]/opts.he)+1
#for unstructured
domain.MeshOptions.he = opts.he
st.assembleDomain(domain)
domain.MeshOptions.triangleOptions = "VApq30Dena%8.8f" % ((opts.he ** 2)/2.0,)

############################################
# ***** Create myTwoPhaseFlowProblem ***** #
############################################
outputStepping = TpFlow.OutputStepping(opts.final_time,dt_output=opts.dt_output,dt_init=0.0001)

theProblem = TpFlow.TwoPhaseFlowProblem(he=opts.he,outputStepping=outputStepping,
                                             domain=domain)

theProblem.addModel(TpFlow.Parameters.ParametersModelRANS2P(),'flow')
theProblem.addModel(TpFlow.Parameters.ParametersModelVOF(),'vof')
theProblem.addModel(TpFlow.Parameters.ParametersModelNCLS(),'ls')
theProblem.addModel(TpFlow.Parameters.ParametersModelRDLS(),'rd')
theProblem.addModel(TpFlow.Parameters.ParametersModelMCorr(),'mcorr')

theProblem.modelDict['flow'].p.initialConditions['p']=zero()
theProblem.modelDict['flow'].p.initialConditions['u']=zero()
theProblem.modelDict['flow'].p.initialConditions['v']=zero()
theProblem.modelDict['vof'].p.initialConditions['vof'] = VF_IC()
theProblem.modelDict['ls'].p.initialConditions['phi'] = PHI_IC()
theProblem.modelDict['rd'].p.initialConditions['phid'] = PHI_IC()
theProblem.modelDict['mcorr'].p.initialConditions['phiCorr'] = zero()
zhang-alvin commented 4 years ago

Test reopening for sizecheck app

zhang-alvin commented 4 years ago

example dambreak case with modified API

(posted above)

zhang-alvin commented 4 years ago

For reference:

git diff upstream/master -- proteus/mprans/SW2DCV.h

yields no changes

ejtovar commented 4 years ago

I know this PR is only for the develop/TPF branch, but why are a lot of the tests failing? @zhang-alvin

cekees commented 4 years ago

I know this PR is only for the develop/TPF branch, but why are a lot of the tests failing? @zhang-alvin

It looks like the main PR changes the TwoPhaseFlowProblem init signature in a way the breaks most of the TwoPhaseFlow tests using the app input. TypeError: __init__() got an unexpected keyword argument 'ns_model'

ejtovar commented 4 years ago

I see, but it's also breaking some of the SWFlow tests.

cekees commented 4 years ago

I see, but it's also breaking some of the SWFlow tests.

It's probably worth your taking a look at the travis output, but you may have to look at the main PR to make sense of what has changed and is actually breaking them. For example that change referenced above doesn't show up in the diffs shown on the PR.

zhang-alvin commented 4 years ago

I know this PR is only for the develop/TPF branch, but why are a lot of the tests failing? @zhang-alvin

It looks like the main PR changes the TwoPhaseFlowProblem init signature in a way the breaks most of the TwoPhaseFlow tests using the app input. TypeError: __init__() got an unexpected keyword argument 'ns_model'

@ejtovar @cekees That's right. The goal here is really to get to an interface that we want and then fix all of the tests at once. An iterative approach would be painfully slow.

I want to also separate the mesh generation from NumericalSolution.py this week for at least TwoPhaseFlow and hopefully streamline some of the physics modules. So if there's no major disagreements with the interface, I'm going to merge this in to develop/TPF tomorrow morning and proceed.