Python proteus.Comm类代码示例

def test_line_gauge_output():
    filename = 'test_line_output.csv'
    silent_rm(filename)

    lines = (((0, 0, 0), (1, 1, 1)),)
    fields = ('u0', )

    l = LineGauges(gauges=((fields, lines),),
                   fileName=filename)

    time_list=[0.0, 1.0, 2.0]
    # good hard-code for up to 1,024 processes or so, but slow (about 8 seconds on single process).
    run_gauge(l, time_list, total_nodes=2048)
    correct_gauge_names = ['u0 [        0         0         0]', 'u0 [ 0.076923  0.076923  0.076923]',
                           'u0 [  0.15385   0.15385   0.15385]', 'u0 [  0.23077   0.23077   0.23077]',
                           'u0 [  0.30769   0.30769   0.30769]', 'u0 [  0.38462   0.38462   0.38462]',
                           'u0 [  0.46154   0.46154   0.46154]', 'u0 [  0.53846   0.53846   0.53846]',
                           'u0 [  0.61538   0.61538   0.61538]', 'u0 [  0.69231   0.69231   0.69231]',
                           'u0 [  0.76923   0.76923   0.76923]', 'u0 [  0.84615   0.84615   0.84615]',
                           'u0 [  0.92308   0.92308   0.92308]', 'u0 [        1         1         1]']
    correct_data = np.asarray([[0., 0., 8.53846154, 17.07692308, 25.61538462,
                                34.15384615, 42.69230769, 51.23076923, 59.76923077, 68.30769231,
                                76.84615385, 85.38461538, 93.92307692, 102.46153846, 111.],
                               [1., 0., 17.07692308, 34.15384615, 51.23076923,
                                68.30769231, 85.38461538, 102.46153846, 119.53846154, 136.61538462,
                                153.69230769, 170.76923077, 187.84615385, 204.92307692, 222.],
                               [2., 0., 25.61538462, 51.23076923, 76.84615385,
                                102.46153846, 128.07692308, 153.69230769, 179.30769231, 204.92307692,
                                230.53846154, 256.15384615, 281.76923077, 307.38461538, 333.]])

    # synchronize processes before attempting to read file

    Comm.get().barrier()

    gauge_names, data = parse_gauge_output(filename)

    eq_(correct_gauge_names, gauge_names)
    npt.assert_allclose(correct_data, data)
    delete_file(filename)

def test_2D_line_integral_gauge_output():
    filename = 'test_2D_line_integral_gauge_output.csv'
    silent_rm(filename)

    lines = (((0, 0, 0), (1, 1, 0)),
             ((0, 0, 0), (1, 0, 0)),
             ((0, 0, 0), (0, 1, 0)),
             ((0, 0.5, 0), (1, 0.5, 0)),
             ((0, 1, 0), (1, 1, 0)),
             ((0.5, 0, 0), (0.5, 1, 0)))
    fields = ('u0', )

    l = LineIntegralGauges(gauges=((fields, lines),),
                   fileName=filename)

    time_list=[0.0, 1.0, 2.0, 2.5]
    run_gauge(l, time_list)

    correct_gauge_names = ['u0 [        0         0         0] - [        1         1         0]',
                           'u0 [        0         0         0] - [        1         0         0]',
                           'u0 [        0         0         0] - [        0         1         0]',
                           'u0 [        0       0.5         0] - [        1       0.5         0]',
                           'u0 [        0         1         0] - [        1         1         0]',
                           'u0 [      0.5         0         0] - [      0.5         1         0]']

    correct_data = np.asarray([[0.,   7.77817459,  0.5,   5.,   5.5,  10.5,   5.5],
                               [1.,  15.55634919,  1.,   10.,  11.,   21.,   11.],
                               [2.,  23.33452378,  1.5,  15.,  16.5,  31.5,  16.5],
                               [2.5, 27.22361108,  1.75, 17.5, 19.25, 36.75, 19.25]])

    # synchronize processes before attempting to read file

    Comm.get().barrier()

    gauge_names, data = parse_gauge_output(filename)
    eq_(correct_gauge_names, gauge_names)

    npt.assert_allclose(correct_data, data)
    delete_file(filename)

def init_mpi_petsc(opts):
    log("Initializing MPI")
    if opts.petscOptions != None:
        petsc_argv = sys.argv[:1]+opts.petscOptions.split()
        log("PETSc options from commandline")
        log(str(petsc_argv))
    else:
        petsc_argv=sys.argv[:1]
    if opts.petscOptionsFile != None:
        petsc_argv=[sys.argv[0]]
        petsc_argv += open(opts.petscOptionsFile).read().split()
        log("PETSc options from commandline")
        log(str(petsc_argv))
    return Comm.init(argv=petsc_argv)

def test_3DparallelLoadPUMI(verbose=0):
    """Test to load 3D parallel PUMI model and mesh"""
    comm = Comm.init()
    eq(comm.size(),2)
    testDir=os.path.dirname(os.path.abspath(__file__))
    domain = Domain.PUMIDomain()
    Model=testDir+ '/Prism.dmg'
    Mesh=testDir + '/Prism.smb'
    domain.PUMIMesh=MeshAdaptPUMI.MeshAdaptPUMI()
    domain.PUMIMesh.loadModelAndMesh(Model, Mesh)
    mesh = MeshTools.TetrahedralMesh()
    mesh.cmesh = cmeshTools.CMesh()
    mesh.convertFromPUMI(domain.PUMIMesh, domain.faceList, domain.regList,parallel = comm.size() > 1, dim = domain.nd)
    eq(mesh.nElements_global,8148)
    eq(mesh.nNodes_global,1880)
    eq(mesh.nEdges_global,11001)
    eq(mesh.nElementBoundaries_global,17270)

def test_2DparallelLoadPUMI(verbose=0):
    """Test to load 2D parallel PUMI model and mesh"""
    comm = Comm.init()
    eq(comm.size(),2)
    testDir=os.path.dirname(os.path.abspath(__file__))
    domain = Domain.PUMIDomain(dim=2)
    Model=testDir+ '/Rectangle.dmg'
    Mesh=testDir + '/Rectangle.smb'
    domain.PUMIMesh=MeshAdaptPUMI.MeshAdaptPUMI()
    domain.PUMIMesh.loadModelAndMesh(Model, Mesh)
    mesh = MeshTools.TriangularMesh()
    mesh.cmesh = cmeshTools.CMesh()
    mesh.convertFromPUMI(domain.PUMIMesh, domain.faceList, domain.regList,parallel = comm.size() > 1, dim = domain.nd)
    eq(mesh.nElements_global,8)
    eq(mesh.nNodes_global,10)
    eq(mesh.nEdges_global,17)
    eq(mesh.nElementBoundaries_global,17)

def test_2DmultiRegion(verbose=0):
    """Test for loading gmsh mesh through PUMI with multiple-regions"""
    testDir=os.path.dirname(os.path.abspath(__file__))
    Model=testDir + '/TwoQuads.dmg'
    Mesh=testDir + '/TwoQuads.smb'
    domain = Domain.PUMIDomain(dim=2) #initialize the domain
    domain.PUMIMesh=MeshAdaptPUMI.MeshAdaptPUMI()
    domain.PUMIMesh.loadModelAndMesh(Model, Mesh)
    domain.faceList=[[14],[12],[11],[13],[15],[16]]
    domain.regList=[[41],[42]]

    mesh = MeshTools.TriangularMesh()
    mesh.cmesh = cmeshTools.CMesh()
    comm = Comm.init()

    mesh.convertFromPUMI(domain.PUMIMesh, domain.faceList,domain.regList, parallel = comm.size() > 1, dim = domain.nd)
    ok(mesh.elementMaterialTypes[0]==1)
    ok(mesh.elementMaterialTypes[-1]==2)

def test_wdot():
    """
    test_wdot

    Verifies that the proteus.LinearAlgebraTools.wdot computes dot product correctly.
    """
    import numpy as np
    import numpy.testing as npt
    from proteus.LinearAlgebraTools import wDot
    from proteus import Comm
    comm = Comm.init()
    x = np.array([-1, 2, 3])
    y = np.array([5, 6, 10])
    h = np.array([0.5, 1.2, 6.0])
    t1 = 191.9        
    t2 = wDot(x, y, h)
    test = npt.assert_almost_equal
    test.description = 'test_wdot'
    assert t1 == approx(t2)

def gauge_setup(nd, total_nodes=None):
    comm = Comm.get()

    #Simplified Physics
    p.name="test_gauges"

    p.nd = nd

    class LinearSolution:
        def uOfXT(self,x,t):
            return (x[0]+10*x[1]+100*x[2])*(t+1.0)

    p.initialConditions = {0:LinearSolution()}
    p.dirichletConditions = {0: lambda x,flag: None}
    p.domain = Domain.RectangularDomain(name="test_gauges_domain")
    p.coefficients = TransportCoefficients.PoissonEquationCoefficients(aOfX = [lambda x: np.eye(p.nd, p.nd)],
                                                                       fOfX = [lambda x: 0], nc=1, nd=p.nd)
    #Simplified and Incomplete Numerics
    n.femSpaces = {0:FemTools.C0_AffineLinearOnSimplexWithNodalBasis}
    n.elementQuadrature = Quadrature.SimplexGaussQuadrature(p.nd,3)
    n.elementBoundaryQuadrature = Quadrature.SimplexGaussQuadrature(p.nd-1,3)
    n.numericalFluxType = NumericalFlux.NoFlux
    n.cfluxtag = None
    n.conservativeFlux = None

    if total_nodes is None:
        total_nodes = 2*comm.size()

    if p.nd == 1:
        mlMesh = build1DMesh(p, total_nodes+1)
    elif p.nd == 2:
        nnx = nny = int(ceil(sqrt(total_nodes)))+1
        mlMesh = build2DMesh(p, nnx, nny)
    elif p.nd == 3:
        nnx = nny = nnz = int(ceil(pow(total_nodes, 1.0/3.0)))+1
        mlMesh = build3DMesh(p, nnx, nny, nnz)

    model = Transport.MultilevelTransport(p, n, mlMesh)

    return model, p.initialConditions

#!/usr/bin/env python
"""
Test module for linear boundary value problems (serial)

This module solves equations of the form

.. math::

  \nabla \cdot \left( a(x) \nabla u \right) = f(x)

"""
from proteus.iproteus import *
from proteus import Comm
comm = Comm.get()
Profiling.logLevel=7
Profiling.verbose=True
import numpy.testing as npt
from nose.tools import ok_ as ok
from nose.tools import eq_ as eq

def test_c0p1():
    import poisson_het_2d_p
    import poisson_het_2d_c0pk_n
    pList = [poisson_het_2d_p]
    nList = [poisson_het_2d_c0pk_n]    
    so = default_so
    so.name = pList[0].name = "poisson_2d_c0p1"+"pe"+`comm.size()`
    so.sList=[default_s]
    opts.logLevel=7
    opts.verbose=True
    opts.profile=True

    petsc_argv=[sys.argv[0]]
    with open(opts.petscOptionsFile) as petsc_file:
        data = petsc_file.readlines()
    def strip_comments(line):
        if '#' in line:
            line = line[:line.index('#')]
        return line

    stripped_data = [strip_comments(line) for line in data]
    petsc_argv += '\n'.join(stripped_data).split()
    log("PETSc options from commandline")
    log(str(petsc_argv))

Comm.argv = petsc_argv
comm = Comm.init()
Profiling.procID=comm.rank()
if opts.logLevel > 0:
    #mwf looks like just gets .py.log if no sso file given
    #Profiling.openLog(args[0][-3:]+".log",opts.logLevel)
    if len(args) == 0:
        Profiling.openLog("proteus.log",opts.logLevel)
    elif len(args[0].split(' ')) == 1:
        Profiling.openLog(args[0].split('.')[0]+".log",opts.logLevel)
    else:
        Profiling.openLog(args[0][-3:]+".log",opts.logLevel)
if opts.logAllProcesses:
    Profiling.logAllProcesses = True

#blanket import statements can go below here now that petsc4py should be initialized
from proteus import *

from proteus import Context
from proteus import Comm
comm = Comm.get()
ctx = Context.get()


# simulation flags for error analysis
#
# simFlagsList is initialized in proteus.iproteus
#

simFlagsList[0]['errorQuantities']=['u']
simFlagsList[0]['errorTypes']= ['numericalSolution'] #compute error in soln and glob. mass bal
simFlagsList[0]['errorNorms']= ['L1','LI'] #compute L2 norm in space or H0 or ...
simFlagsList[0]['errorTimes']= ['Last'] #'All', 'Last'
simFlagsList[0]['echo']=True
simFlagsList[0]['echoRelativeErrors']=True
#
start
quit

    def __init__(self,
                 uDict,
                 phiDict,
                 testSpaceDict,
                 matType,
                 dofBoundaryConditionsDict,
                 dofBoundaryConditionsSetterDict,
                 coefficients,
                 elementQuadrature,
                 elementBoundaryQuadrature,
                 fluxBoundaryConditionsDict=None,
                 advectiveFluxBoundaryConditionsSetterDict=None,
                 diffusiveFluxBoundaryConditionsSetterDictDict=None,
                 stressTraceBoundaryConditionsSetterDict=None,
                 stabilization=None,
                 shockCapturing=None,
                 conservativeFluxDict=None,
                 numericalFluxType=None,
                 TimeIntegrationClass=None,
                 massLumping=False,
                 reactionLumping=False,
                 options=None,
                 name='defaultName',
                 reuse_trial_and_test_quadrature=True,
                 sd=True,
                 movingDomain=False,
                 bdyNullSpace=True):
        self.bdyNullSpace = bdyNullSpace
        from proteus import Comm
        #
        # set the objects describing the method and boundary conditions
        #
        self.movingDomain = movingDomain
        self.tLast_mesh = None
        #
        self.name = name
        self.sd = sd
        self.Hess = False
        self.lowmem = True
        self.timeTerm = True  # allow turning off  the  time derivative
        # self.lowmem=False
        self.testIsTrial = True
        self.phiTrialIsTrial = True
        self.u = uDict
        self.ua = {}  # analytical solutions
        self.phi = phiDict
        self.dphi = {}
        self.matType = matType
        # mwf try to reuse test and trial information across components if
        # spaces are the same
        self.reuse_test_trial_quadrature = reuse_trial_and_test_quadrature  # True#False
        if self.reuse_test_trial_quadrature:
            for ci in range(1, coefficients.nc):
                assert self.u[ci].femSpace.__class__.__name__ == self.u[
                    0].femSpace.__class__.__name__, "to reuse_test_trial_quad all femSpaces must be the same!"
        # Simplicial Mesh
        # assume the same mesh for  all components for now
        self.mesh = self.u[0].femSpace.mesh
        self.testSpace = testSpaceDict
        self.dirichletConditions = dofBoundaryConditionsDict
        # explicit Dirichlet  conditions for now, no Dirichlet BC constraints
        self.dirichletNodeSetList = None
        self.coefficients = coefficients
        self.coefficients.initializeMesh(self.mesh)
        self.nc = self.coefficients.nc
        self.stabilization = stabilization
        self.shockCapturing = shockCapturing
        # no velocity post-processing for now
        self.conservativeFlux = conservativeFluxDict
        self.fluxBoundaryConditions = fluxBoundaryConditionsDict
        self.diffusiveFluxBoundaryConditionsSetterDictDict = diffusiveFluxBoundaryConditionsSetterDictDict
        # determine whether  the stabilization term is nonlinear
        self.stabilizationIsNonlinear = False
        # cek come back
        if self.stabilization is not None:
            for ci in range(self.nc):
                if coefficients.mass.has_key(ci):
                    for flag in coefficients.mass[ci].values():
                        if flag == 'nonlinear':
                            self.stabilizationIsNonlinear = True
                if coefficients.advection.has_key(ci):
                    for flag in coefficients.advection[ci].values():
                        if flag == 'nonlinear':
                            self.stabilizationIsNonlinear = True
                if coefficients.diffusion.has_key(ci):
                    for diffusionDict in coefficients.diffusion[ci].values():
                        for flag in diffusionDict.values():
                            if flag != 'constant':
                                self.stabilizationIsNonlinear = True
                if coefficients.potential.has_key(ci):
                    for flag in coefficients.potential[ci].values():
                        if flag == 'nonlinear':
                            self.stabilizationIsNonlinear = True
                if coefficients.reaction.has_key(ci):
                    for flag in coefficients.reaction[ci].values():
                        if flag == 'nonlinear':
                            self.stabilizationIsNonlinear = True
                if coefficients.hamiltonian.has_key(ci):
                    for flag in coefficients.hamiltonian[ci].values():
                        if flag == 'nonlinear':
#.........这里部分代码省略.........

log = Profiling.logEvent

log("Initializing MPI")#this won't show up unless logAllProcesses=True because we don't know our procID yet
if opts.petscOptions != None:
    petsc_argv = sys.argv[:1]+opts.petscOptions.split()
    log("PETSc options from commandline")
    log(str(petsc_argv))
else:
    petsc_argv=sys.argv[:1]
if opts.petscOptionsFile != None:
    petsc_argv=[sys.argv[0]]
    petsc_argv += open(opts.petscOptionsFile).read().split()
    log("PETSc options from commandline")
    log(str(petsc_argv))

comm = Comm.init(argv=petsc_argv)
#blanket import statements can go below here now that petsc4py should be initialized
Profiling.procID=comm.rank()
log("Initialized MPI")
if opts.viewer is not False:
    from proteusGraphical import *
from proteus import *

log("Adding "+str(opts.probDir)+" to path for loading modules")
probDir = str(opts.probDir)
if probDir not in sys.path:
    sys.path.insert(0,probDir)

pList = []
nList = []
sList = []

def test_2DgmshLoadAndAdapt(verbose=0):
    """Test for loading gmsh mesh through PUMI, estimating error and adapting for 
    a 2D Couette flow case"""
    testDir=os.path.dirname(os.path.abspath(__file__))
    Model=testDir + '/Couette2D.null'
    Mesh=testDir + '/Couette2D.msh'
    domain = Domain.PUMIDomain(dim=2) #initialize the domain
    domain.PUMIMesh=MeshAdaptPUMI.MeshAdaptPUMI(hmax=0.01, hmin=0.008, numIter=1,sfConfig='ERM',maType='isotropic',targetError=1)
    domain.PUMIMesh.loadModelAndMesh(Model, Mesh)
    domain.faceList=[[14],[12],[11],[13]]

    mesh = MeshTools.TriangularMesh()
    mesh.cmesh = cmeshTools.CMesh()
    comm = Comm.init()

    nElements_initial = mesh.nElements_global
    mesh.convertFromPUMI(domain.PUMIMesh, domain.faceList,domain.regList, parallel = comm.size() > 1, dim = domain.nd)

    domain.PUMIMesh.transferFieldToPUMI("coordinates",mesh.nodeArray)

    rho = numpy.array([998.2,998.2])
    nu = numpy.array([1.004e-6, 1.004e-6])
    g = numpy.asarray([0.0,0.0])
    deltaT = 1.0 #dummy number
    domain.PUMIMesh.transferPropertiesToPUMI(rho,nu,g,deltaT)

    #Couette Flow
    Lz = 0.05
    Uinf = 2e-3
    #hard code solution
    vector=numpy.zeros((mesh.nNodes_global,3),'d')
    dummy = numpy.zeros(mesh.nNodes_global); 
    vector[:,0] = Uinf*mesh.nodeArray[:,1]/Lz #v-velocity
    vector[:,1] = dummy
    vector[:,2] = dummy
    domain.PUMIMesh.transferFieldToPUMI("velocity", vector)
    del vector
    del dummy

    scalar=numpy.zeros((mesh.nNodes_global,1),'d')
    domain.PUMIMesh.transferFieldToPUMI("p", scalar)

    scalar[:,0] = mesh.nodeArray[:,1]
    domain.PUMIMesh.transferFieldToPUMI("phi", scalar)
    del scalar

    scalar = numpy.zeros((mesh.nNodes_global,1),'d')+1.0
    domain.PUMIMesh.transferFieldToPUMI("vof", scalar)

    errorTotal=domain.PUMIMesh.get_local_error()
    ok(errorTotal<1e-14)

    ok(domain.PUMIMesh.willAdapt(),1)

    domain.PUMIMesh.adaptPUMIMesh()
    
    mesh = MeshTools.TriangularMesh()
    mesh.convertFromPUMI(domain.PUMIMesh,
                     domain.faceList,
                     domain.regList,
                     parallel = comm.size() > 1,
                     dim = domain.nd)
    nElements_final = mesh.nElements_global
    ok(nElements_final>nElements_initial)

#import sys
#from proteus import Comm
import os
import sys
import socket
import cPickle
import numpy
import proteus
#remove blankent import statements until after Comm initialized for petsc4py
# from proteus import *
from proteus import Profiling,Comm
from warnings import *
import optparse
import sys
import pstats
import pdb
import petsc4py
Profiling.openLog("proteus.log",7)
print sys.argv[:1]
comm = Comm.init(argv=sys.argv[:1])
comm = Comm.init(argv=sys.argv[:1])
#petsc4py.init(sys.argv[:1])
print comm.rank(),comm.size()
print "Hellow World from",comm.rank()
from proteus import *
Profiling.procID=comm.rank()