1 | import numpy as np
|
---|
2 | import os
|
---|
3 | from ContourToMesh import ContourToMesh
|
---|
4 | import MatlabFuncs as m
|
---|
5 | import PythonFuncs as p
|
---|
6 |
|
---|
7 |
|
---|
8 | def FlagElements(md, region):
|
---|
9 | """
|
---|
10 | FLAGELEMENTS - flag the elements in an region
|
---|
11 |
|
---|
12 | The region can be given with an exp file, a list of elements or vertices
|
---|
13 |
|
---|
14 | Usage:
|
---|
15 | flag = FlagElements(md, region)
|
---|
16 |
|
---|
17 | Example:
|
---|
18 | flag = FlagElements(md, 'all')
|
---|
19 | flag = FlagElements(md, '')
|
---|
20 | flag = FlagElements(md, 'Domain.exp')
|
---|
21 | flag = FlagElements(md, '~Domain.exp')
|
---|
22 | """
|
---|
23 |
|
---|
24 | if isinstance(region, str):
|
---|
25 | if not region:
|
---|
26 | flag = np.zeros(md.mesh.numberofelements, bool)
|
---|
27 | invert = 0
|
---|
28 | elif m.strcmpi(region, 'all'):
|
---|
29 | flag = np.ones(md.mesh.numberofelements, bool)
|
---|
30 | invert = 0
|
---|
31 | else:
|
---|
32 | #make sure that we actually don't want the elements outside the domain outline!
|
---|
33 | if m.strcmpi(region[0], '~'):
|
---|
34 | region = region[1:]
|
---|
35 | invert = 1
|
---|
36 | else:
|
---|
37 | invert = 0
|
---|
38 |
|
---|
39 | #does the region domain outline exist or do we have to look for xlim, ylim in basinzoom?
|
---|
40 | if not os.path.exists(region):
|
---|
41 | if len(region) > 3 and not m.strcmp(region[-4:], '.exp'):
|
---|
42 | raise IOError("Error: File 'region' not found!" % region)
|
---|
43 | raise RuntimeError("FlagElements.py calling basinzoom.py is not complete.")
|
---|
44 | xlim, ylim = basinzoom('basin', region)
|
---|
45 | flag_nodes = p.logical_and_n(md.mesh.x < xlim[1], md.mesh.x > xlim[0], md.mesh.y < ylim[1], md.mesh.y > ylim[0])
|
---|
46 | flag = np.prod(flag_nodes[md.mesh.elements], axis=1).astype(bool)
|
---|
47 | else:
|
---|
48 | #ok, flag elements
|
---|
49 | flag = ContourToMesh(md.mesh.elements[:, 0:3].copy(), md.mesh.x, md.mesh.y, region, 'element', 1)
|
---|
50 | flag = flag.astype(bool)
|
---|
51 |
|
---|
52 | if invert:
|
---|
53 | flag = np.logical_not(flag)
|
---|
54 |
|
---|
55 | elif isinstance(region, np.ndarray) or isinstance(region, bool):
|
---|
56 | if np.size(region, 0) == md.mesh.numberofelements:
|
---|
57 | flag = region
|
---|
58 | elif np.size(region, 0) == md.mesh.numberofvertices:
|
---|
59 | flag = (np.sum(region[md.mesh.elements - 1] > 0, axis=1) == np.size(md.mesh.elements, 1))
|
---|
60 | else:
|
---|
61 | raise TypeError("Flaglist for region must be of same size as number of elements in model.")
|
---|
62 |
|
---|
63 | else:
|
---|
64 | raise TypeError("Invalid region option")
|
---|
65 |
|
---|
66 | return flag
|
---|