Changeset 13520


Ignore:
Timestamp:
10/03/12 10:38:28 (12 years ago)
Author:
jschierm
Message:

FIX: Fix handling of variable arguments (and more).

File:
1 edited

Legend:

Unmodified
Added
Removed
  • issm/trunk-jpl/src/m/mesh/triangle.py

    r13009 r13520  
    1 from numpy import *
    2 import TriMesh as tm
    3 import NodeConnectivity as nc
    4 import ElementConnectivity as ec
     1import numpy
     2from TriMesh import *
     3from NodeConnectivity import *
     4from ElementConnectivity import *
     5from MatlabFuncs import *
    56
    6 def triangle(md, domainname, resolution,riftname=''):
    7         #TRIANGLE - create model mesh using the triangle package
    8         #
    9         #   This routine creates a model mesh using TriMesh and a domain outline, to within a certain resolution
    10         #   where md is a @model object, domainname is the name of an Argus domain outline file,
    11         #   and resolution is a characteristic length for the mesh (same unit as the domain outline
    12         #   unit). Riftname is an optional argument (Argus domain outline) describing rifts.
    13         #
    14         #   Usage:
    15         #      md=triangle(md,domainname,resolution)
    16         #   or md=triangle(md,domainname, resolution, riftname)
    17         #
    18         #   Examples:
    19         #      md=triangle(md,'DomainOutline.exp',1000);
    20         #      md=triangle(md,'DomainOutline.exp',1000,'Rifts.exp');
     7def triangle(md,domainname,*args):
     8        """
     9        TRIANGLE - create model mesh using the triangle package
    2110
     11           This routine creates a model mesh using TriMesh and a domain outline, to within a certain resolution
     12           where md is a @model object, domainname is the name of an Argus domain outline file,
     13           and resolution is a characteristic length for the mesh (same unit as the domain outline
     14           unit). Riftname is an optional argument (Argus domain outline) describing rifts.
     15
     16           Usage:
     17              md=triangle(md,domainname,resolution)
     18           or md=triangle(md,domainname, resolution, riftname)
     19
     20           Examples:
     21              md=triangle(md,'DomainOutline.exp',1000);
     22              md=triangle(md,'DomainOutline.exp',1000,'Rifts.exp');
     23        """
    2224
    2325        #Figure out a characteristic area. Resolution is a node oriented concept (ex a 1000m  resolution node would
    2426        #be made of 1000*1000 area squares).
    2527
     28        if len(args)==1:
     29                resolution=args[0]
     30                riftname=''
     31        if len(args)==2:
     32                riftname=args[0]
     33                resolution=args[1]
     34
    2635        #Check that mesh was not already run, and warn user:
    27         if md.mesh.numberofelements != 0.:
    28                 choice = input('This model already has a mesh. Are you sure you want to go ahead? (y/n)')
    29                 if choice != 'y':
     36        if md.mesh.numberofelements:
     37                choice = raw_input('This model already has a mesh. Are you sure you want to go ahead? (y/n)')
     38                if not strcmp(choice,'y'):
    3039                        print 'no meshing done ... exiting'
    31                         return []
    32                
    33         area = resolution**2.
     40                        return None
     41
     42        area = resolution**2
    3443
    3544        #Mesh using TriMesh
    36         [md.mesh.elements,md.mesh.x,md.mesh.y,md.mesh.segments,md.mesh.segmentmarkers]=tm.TriMesh(domainname,riftname,area)
    37 
     45        [md.mesh.elements,md.mesh.x,md.mesh.y,md.mesh.segments,md.mesh.segmentmarkers]=TriMesh(domainname,riftname,area)
    3846
    3947        #Fill in rest of fields:
    40         md.mesh.numberofelements = len(md.mesh.elements)
    41         md.mesh.numberofvertices = len(md.mesh.x)
    42         md.mesh.z = zeros(md.mesh.numberofvertices)
    43         md.mesh.vertexonboundary = zeros(md.mesh.numberofvertices)
     48        md.mesh.numberofelements = numpy.size(md.mesh.elements,axis=0)
     49        md.mesh.numberofvertices = numpy.size(md.mesh.x)
     50        md.mesh.z = numpy.zeros(md.mesh.numberofvertices)
     51        md.mesh.vertexonboundary = numpy.zeros(md.mesh.numberofvertices)
    4452        md.mesh.vertexonboundary[md.mesh.segments[:,0:2].astype(int)-1] = 1.
    45         md.mesh.vertexonbed = ones(md.mesh.numberofvertices)
    46         md.mesh.vertexonsurface = ones(md.mesh.numberofvertices)
    47         md.mesh.elementonbed = ones(md.mesh.numberofelements)
    48         md.mesh.elementonsurface = ones(md.mesh.numberofelements)
     53        md.mesh.vertexonbed = numpy.ones(md.mesh.numberofvertices)
     54        md.mesh.vertexonsurface = numpy.ones(md.mesh.numberofvertices)
     55        md.mesh.elementonbed = numpy.ones(md.mesh.numberofelements)
     56        md.mesh.elementonsurface = numpy.ones(md.mesh.numberofelements)
    4957
    5058        #Now, build the connectivity tables for this mesh.
    51         [md.mesh.vertexconnectivity]= nc.NodeConnectivity(md.mesh.elements, md.mesh.numberofvertices)
    52         [md.mesh.elementconnectivity] = ec.ElementConnectivity(md.mesh.elements, md.mesh.vertexconnectivity)
    53        
     59        [md.mesh.vertexconnectivity]= NodeConnectivity(md.mesh.elements, md.mesh.numberofvertices)
     60        [md.mesh.elementconnectivity] = ElementConnectivity(md.mesh.elements, md.mesh.vertexconnectivity)
     61
    5462        #type of model
    5563        md.mesh.dimension = 2.
     64
    5665        return md
     66
Note: See TracChangeset for help on using the changeset viewer.