Secure your code as it's written. Use Snyk Code to scan source code in minutes - no build needed - and fix issues immediately.
# Set size and position of dense sphere.
sphereRadius = 0.1
sphereCentre = (0., 0.7)
# Create mesh and finite element variables
# ------
# In[3]:
mesh = uw.mesh.FeMesh_Cartesian( elementType = ("Q1/dQ0"),
elementRes = (res, res),
minCoord = (-1., 0.),
maxCoord = (1., 1.))
velocityField = uw.mesh.MeshVariable( mesh=mesh, nodeDofCount=dim )
pressureField = uw.mesh.MeshVariable( mesh=mesh.subMesh, nodeDofCount=1 )
# Set initial conditions and boundary conditions
# ----------
#
# **Initial and boundary conditions**
#
# Initialise the velocity and pressure fields to zero.
# In[4]:
velocityField.data[:] = [0.,0.]
pressureField.data[:] = 0.
import underworld as uw
import math
from underworld import function as fn
linearMesh = uw.mesh.FeMesh_Cartesian("Q2/dQ1", (32,32), (0.,0.), (1.,1.))
constantMesh = linearMesh.subMesh
velocityField = uw.mesh.MeshVariable(linearMesh,2)
velocityField.data[:] = (0.,0.)
pressureField = uw.mesh.MeshVariable(constantMesh,1)
pressureField.data[:] = 0.
# freeslip
IWalls = linearMesh.specialSets["MinI_VertexSet"] + linearMesh.specialSets["MaxI_VertexSet"]
JWalls = linearMesh.specialSets["MinJ_VertexSet"] + linearMesh.specialSets["MaxJ_VertexSet"]
freeslip = uw.conditions.DirichletCondition(velocityField, (IWalls, JWalls))
solA = fn.analytic.SolCx()
stokesSystem = uw.systems.Stokes(velocityField,pressureField,solA.fn_viscosity,solA.fn_bodyforce,conditions=[freeslip,], rtolerance=1.e-5)
# To set up Uzawa solver with mumps
def __init__(self, velocityMeshVariable, pressureMeshVariable, temperatureMeshVariable=None, pic=True, bodyforces=[], **kwargs):
# build parent
super(MeshVariable,self).__init__(**kwargs)
# now create a bunch of stuff
Ra=1.0
if temperatureMeshVariable is None:
#then build a temperature field
mesh = velocityMeshVariable._mesh
temperatureMeshVariable = uw.mesh.MeshVariable( mesh=mesh, nodeDofCount=1, dataType="double")
#create our own int swarm here?
thermalBuoyancyTerm = ThermalBuoyancy(temperatureMeshVariable, intswarm, Ra)
f = ForceVector(velocityMeshVariable,[])
h = ForceVector(pressureMeshVariable,[])
u = SolutionVector(velocityMeshVariable)
p = SolutionVector(pressureMeshVariable)
# make directories if they don't exist
if not os.path.isdir(outputPath):
os.makedirs(outputPath)
if not os.path.isdir(filePath):
os.makedirs(filePath)
# Create mesh and finite element variables
# ------
#
# Set up mesh and field variables that are solved on the mesh. See user guides for details of this process.
# In[ ]:
mesh = uw.mesh.FeMesh_Cartesian( elementType = ("Q1/dQ0"),
elementRes = (res, res),
minCoord = (0., 0.),
maxCoord = (1., 1.))
velocityField = uw.mesh.MeshVariable( mesh=mesh, nodeDofCount=2 )
pressureField = uw.mesh.MeshVariable( mesh=mesh.subMesh, nodeDofCount=1 )
temperatureField = uw.mesh.MeshVariable( mesh=mesh, nodeDofCount=1 )
temperatureDotField = uw.mesh.MeshVariable( mesh=mesh, nodeDofCount=1 )
# Set initial conditions and boundary conditions
# ----------
#
# **Initial and boundary conditions**
#
# Either set by perturbation function or load data from file.
def __init__(self, meshVariable, **kwargs):
if not isinstance(meshVariable, mesh.MeshVariable):
raise TypeError("'meshVariable' object passed in must be of type 'MeshVariable'")
self._meshVariable = meshVariable
# build parent
super(Vector,self).__init__(**kwargs)
# In[ ]:
LoadFromFile = True
savedRes = 128
# **If loading from file**
#
# Read (``savedRes`` $\times$ ``savedRes``) resolution data for $P$, $v$ and $T$ fields as well as existing summary statistics data. These are converted into lists so that the main time loop below will append with new values.
#
# In[ ]:
if(LoadFromFile == True):
# set up mesh for savedRes*savedRes data file
meshSaved = uw.mesh.FeMesh_Cartesian( elementType = ("Q1/dQ0"),
elementRes = (savedRes, savedRes),
minCoord = (0., 0.),
maxCoord = (1., 1.) )
temperatureFieldSaved = uw.mesh.MeshVariable( mesh=meshSaved, nodeDofCount=1 )
temperatureDotFieldSaved = uw.mesh.MeshVariable( mesh=meshSaved, nodeDofCount=1 )
pressureFieldSaved = uw.mesh.MeshVariable( mesh=meshSaved.subMesh, nodeDofCount=1 )
velocityFieldSaved = uw.mesh.MeshVariable( mesh=meshSaved, nodeDofCount=2 )
temperatureFieldSaved.load( inputPath+'temperatureField.h5' )
temperatureDotFieldSaved.load( inputPath+'temperatureDotField.h5' )
pressureFieldSaved.load( inputPath+'pressureField.h5')
velocityFieldSaved.load( inputPath+'velocityField.h5')
if(res==savedRes): # no remeshing needed, copy data directly
temperatureField.data[:] = temperatureFieldSaved.data[:]
pressureField.data[:] = pressureFieldSaved.data[:]
self.specialSets["Empty"] = lambda selfobject: uw.mesh.FeMesh_IndexSet( object = selfobject,
topologicalIndex = 0,
size = libUnderworld.StgDomain.Mesh_GetDomainSize( selfobject._mesh, libUnderworld.StgDomain.MT_VERTEX ))
outputPath = os.path.join(workdir,"outputs/")
if uw.rank() == 0:
if not os.path.exists(outputPath):
os.makedirs(outputPath)
# Create mesh and finite element variables
# ------
#
# Note: the use of a pressure-sensitive rheology suggests that it is important to use a Q2/dQ1 element
# In[5]:
mesh = uw.mesh.FeMesh_Cartesian( elementType = ("Q1/dQ0"),
elementRes = ( resX, resY),
minCoord = ( xmin, 0.),
maxCoord = ( xmax, 1.),
periodic = [False, False] )
velocityField = uw.mesh.MeshVariable( mesh=mesh, nodeDofCount=mesh.dim )
pressureField = uw.mesh.MeshVariable( mesh=mesh.subMesh, nodeDofCount=1 )
velocityField.data[:] = [0.,0.]
pressureField.data[:] = 0.
# ### Boundary conditions
#
# Pure shear with moving walls — all boundaries are zero traction with
#
# The mesh object has both a primary and sub mesh. "Q1/dQ0" produces a primary mesh with element type Q1 and a sub-mesh with elements type dQ0. Q1 elements have node points at the element corners, dQ0 elements have a single node at the elements centre.
mesh = uw.mesh.FeMesh_Cartesian( elementType = ("Q1/dQ0"),
elementRes = (int(boxLength*res), res),
minCoord = (0., 0.),
maxCoord = (boxLength, boxHeight))
Tmesh = mesh
# Create mesh variables. Note the pressure field uses the sub-mesh, the temperature uses the same mesh at this stage. The `temperatureDotField` variable is a work array for the advection-diffusion solver.
velocityField = uw.mesh.MeshVariable( mesh=mesh, nodeDofCount=dim )
pressureField = uw.mesh.MeshVariable( mesh=mesh.subMesh, nodeDofCount=1 )
temperatureField = uw.mesh.MeshVariable( mesh=Tmesh, nodeDofCount=1 )
temperatureDotField = uw.mesh.MeshVariable( mesh=Tmesh, nodeDofCount=1 )
# Initialise values
velocityField.data[:] = [0.,0.]
pressureField.data[:] = 0.
temperatureDotField.data[:] = 0.
# Set up material parameters and functions
# -----
#
# Set functions for viscosity, density and buoyancy force. These functions and variables only need to be defined at the beginning of the simulation, not each timestep.
nodeDofCount : int
Number of degrees of freedom per node the variable will have
Returns
-------
underworld.mesh.MeshVariable
The newly created mesh variable.
Example
-------
>>> linearMesh = uw.mesh.FeMesh_Cartesian( elementType='Q1/dQ0', elementRes=(16,16), minCoord=(0.,0.), maxCoord=(1.,1.) )
>>> scalarFeVar = linearMesh.add_variable( nodeDofCount=1, dataType="double" )
>>> q0field = linearMesh.subMesh.add_variable( 1 ) # adds variable to secondary elementType discretisation
"""
var = uw.mesh.MeshVariable(self, nodeDofCount, dataType, **kwargs)
return var