[MOAB-dev] commit/MOAB: tautges: Merge branch 'master' into tautges/par_spatial_locator
commits-noreply at bitbucket.org
commits-noreply at bitbucket.org
Thu Feb 13 11:25:21 CST 2014
1 new commit in MOAB:
https://bitbucket.org/fathomteam/moab/commits/fe8ac7589ddd/
Changeset: fe8ac7589ddd
Branch: tautges/par_spatial_locator
User: tautges
Date: 2014-02-13 18:25:03
Summary: Merge branch 'master' into tautges/par_spatial_locator
Affected #: 39 files
diff --git a/.gitignore b/.gitignore
index 482fcd4..e60d5d8 100644
--- a/.gitignore
+++ b/.gitignore
@@ -178,7 +178,6 @@ test/obb/obb_time
test/obb/obb_tree_tool
test/obb_test
test/oldinc/test_oldinc
-test/read_mpas_nc
test/parallel/*.h5m
test/parallel/mbparallelcomm_test
test/parallel/mhdf_parallel
diff --git a/MeshFiles/unittest/Makefile.am b/MeshFiles/unittest/Makefile.am
index a1784ee..ab3e23e 100644
--- a/MeshFiles/unittest/Makefile.am
+++ b/MeshFiles/unittest/Makefile.am
@@ -33,5 +33,6 @@ EXTRA_DIST = 125hex.g \
poly14.txt \
BedCrop2.h5m \
mpas_p8.h5m \
- Homme_2pt.h5m
+ Homme_2pt.h5m \
+ surfrandomtris-4part.h5m
diff --git a/MeshFiles/unittest/io/cube.sat b/MeshFiles/unittest/io/cube.sat
new file mode 100644
index 0000000..58fd508
--- /dev/null
+++ b/MeshFiles/unittest/io/cube.sat
@@ -0,0 +1,117 @@
+1900 0 1 0
+10 Cubit 12.2 17 ACIS 19.0.2 Linux 24 Mon Jan 6 14:16:03 2014
+1 9.9999999999999995e-07 1e-10
+body $1 -1 -1 $-1 $2 $-1 $-1 T -5 -5 -5 5 5 5 #
+integer_attrib-name_attrib-gen-attrib $-1 -1 $-1 $-1 $0 2 1 1 1 1 1 1 1 1 1 1 1 1 1 0 1 1 1 @8 CUBIT_ID 1 #
+lump $3 -1 -1 $-1 $-1 $4 $0 T -5 -5 -5 5 5 5 #
+integer_attrib-name_attrib-gen-attrib $-1 -1 $-1 $-1 $2 2 1 1 1 1 1 1 1 1 1 1 1 1 1 0 1 1 1 @8 CUBIT_ID 1 #
+shell $-1 -1 -1 $-1 $-1 $-1 $5 $-1 $2 T -5 -5 -5 5 5 5 #
+face $6 -1 -1 $-1 $7 $8 $4 $-1 $9 forward single T -5 -5 5 5 5 5 F #
+integer_attrib-name_attrib-gen-attrib $-1 -1 $-1 $-1 $5 2 1 1 1 1 1 1 1 1 1 1 1 1 1 0 1 1 1 @8 CUBIT_ID 1 #
+face $10 -1 -1 $-1 $11 $12 $4 $-1 $13 reversed single T -5 -5 -5 5 5 -5 F #
+loop $-1 -1 -1 $-1 $-1 $14 $5 T -5 -5 5 5 5 5 unknown #
+plane-surface $-1 -1 -1 $-1 0 0 5 0 0 1 1 0 0 forward_v I I I I #
+integer_attrib-name_attrib-gen-attrib $-1 -1 $-1 $-1 $7 2 1 1 1 1 1 1 1 1 1 1 1 1 1 0 1 1 1 @8 CUBIT_ID 2 #
+face $15 -1 -1 $-1 $16 $17 $4 $-1 $18 reversed single T -5 -5 -5 5 -5 5 F #
+loop $-1 -1 -1 $-1 $-1 $19 $7 T -5 -5 -5 5 5 -5 unknown #
+plane-surface $-1 -1 -1 $-1 0 0 -5 0 0 1 1 0 0 forward_v I I I I #
+coedge $-1 -1 -1 $-1 $20 $21 $22 $23 forward $8 $-1 #
+integer_attrib-name_attrib-gen-attrib $-1 -1 $-1 $-1 $11 2 1 1 1 1 1 1 1 1 1 1 1 1 1 0 1 1 1 @8 CUBIT_ID 3 #
+face $24 -1 -1 $-1 $25 $26 $4 $-1 $27 reversed single T -5 -5 -5 -5 5 5 F #
+loop $-1 -1 -1 $-1 $-1 $28 $11 T -5 -5 -5 5 -5 5 unknown #
+plane-surface $-1 -1 -1 $-1 0 -5 0 0 1 -0 -0 0 1 forward_v I I I I #
+coedge $-1 -1 -1 $-1 $29 $30 $31 $32 forward $12 $-1 #
+coedge $-1 -1 -1 $-1 $33 $14 $34 $35 forward $8 $-1 #
+coedge $-1 -1 -1 $-1 $14 $33 $36 $37 forward $8 $-1 #
+coedge $-1 -1 -1 $-1 $38 $39 $14 $23 reversed $40 $-1 #
+edge $41 -1 -1 $-1 $42 -5 $43 5 $22 $44 forward @7 unknown T 5 -5 5 5 5 5 #
+integer_attrib-name_attrib-gen-attrib $-1 -1 $-1 $-1 $16 2 1 1 1 1 1 1 1 1 1 1 1 1 1 0 1 1 1 @8 CUBIT_ID 4 #
+face $45 -1 -1 $-1 $46 $47 $4 $-1 $48 reversed single T -5 5 -5 5 5 5 F #
+loop $-1 -1 -1 $-1 $-1 $49 $16 T -5 -5 -5 -5 5 5 unknown #
+plane-surface $-1 -1 -1 $-1 -5 0 0 1 0 0 0 0 -1 forward_v I I I I #
+coedge $-1 -1 -1 $-1 $50 $36 $51 $52 forward $17 $-1 #
+coedge $-1 -1 -1 $-1 $53 $19 $50 $54 forward $12 $-1 #
+coedge $-1 -1 -1 $-1 $19 $53 $55 $56 forward $12 $-1 #
+coedge $-1 -1 -1 $-1 $39 $38 $19 $32 reversed $40 $-1 #
+edge $57 -1 -1 $-1 $58 -5 $59 5 $31 $60 forward @7 unknown T 5 -5 -5 5 5 -5 #
+coedge $-1 -1 -1 $-1 $21 $20 $61 $62 forward $8 $-1 #
+coedge $-1 -1 -1 $-1 $63 $64 $20 $35 reversed $47 $-1 #
+edge $65 -1 -1 $-1 $43 -5 $66 5 $34 $67 forward @7 unknown T -5 5 5 5 5 5 #
+coedge $-1 -1 -1 $-1 $28 $68 $21 $37 reversed $17 $-1 #
+edge $69 -1 -1 $-1 $70 -5 $42 5 $36 $71 forward @7 unknown T -5 -5 5 5 -5 5 #
+coedge $-1 -1 -1 $-1 $31 $22 $68 $72 forward $40 $-1 #
+coedge $-1 -1 -1 $-1 $22 $31 $63 $73 reversed $40 $-1 #
+loop $-1 -1 -1 $-1 $-1 $38 $46 T 5 -5 -5 5 5 5 unknown #
+integer_attrib-name_attrib-gen-attrib $-1 -1 $-1 $-1 $23 2 1 1 1 1 1 1 1 1 1 1 1 1 1 0 1 1 1 @8 CUBIT_ID 1 #
+vertex $74 -1 -1 $-1 $23 $75 #
+vertex $76 -1 -1 $-1 $23 $77 #
+straight-curve $-1 -1 -1 $-1 5 0 5 0 1 0 I I #
+integer_attrib-name_attrib-gen-attrib $-1 -1 $-1 $-1 $25 2 1 1 1 1 1 1 1 1 1 1 1 1 1 0 1 1 1 @8 CUBIT_ID 5 #
+face $78 -1 -1 $-1 $-1 $40 $4 $-1 $79 reversed single T 5 -5 -5 5 5 5 F #
+loop $-1 -1 -1 $-1 $-1 $63 $25 T -5 5 -5 5 5 5 unknown #
+plane-surface $-1 -1 -1 $-1 0 5 0 0 -1 0 0 0 -1 forward_v I I I I #
+coedge $-1 -1 -1 $-1 $80 $61 $64 $81 forward $26 $-1 #
+coedge $-1 -1 -1 $-1 $68 $28 $29 $54 reversed $17 $-1 #
+coedge $-1 -1 -1 $-1 $61 $80 $28 $52 reversed $26 $-1 #
+edge $82 -1 -1 $-1 $70 -5 $83 5 $51 $84 forward @7 unknown T -5 -5 -5 -5 -5 5 #
+coedge $-1 -1 -1 $-1 $30 $29 $80 $85 forward $12 $-1 #
+edge $86 -1 -1 $-1 $59 -5 $83 5 $50 $87 forward @7 unknown T -5 -5 -5 5 -5 -5 #
+coedge $-1 -1 -1 $-1 $64 $63 $30 $56 reversed $47 $-1 #
+edge $88 -1 -1 $-1 $89 -5 $58 5 $55 $90 forward @7 unknown T -5 5 -5 5 5 -5 #
+integer_attrib-name_attrib-gen-attrib $-1 -1 $-1 $-1 $32 2 1 1 1 1 1 1 1 1 1 1 1 1 1 0 1 1 1 @8 CUBIT_ID 5 #
+vertex $91 -1 -1 $-1 $32 $92 #
+vertex $93 -1 -1 $-1 $72 $94 #
+straight-curve $-1 -1 -1 $-1 5 0 -5 0 -1 0 I I #
+coedge $-1 -1 -1 $-1 $49 $51 $33 $62 reversed $26 $-1 #
+edge $95 -1 -1 $-1 $66 -5 $70 5 $61 $96 forward @7 unknown T -5 -5 5 -5 5 5 #
+coedge $-1 -1 -1 $-1 $55 $34 $39 $73 forward $47 $-1 #
+coedge $-1 -1 -1 $-1 $34 $55 $49 $81 reversed $47 $-1 #
+integer_attrib-name_attrib-gen-attrib $-1 -1 $-1 $-1 $35 2 1 1 1 1 1 1 1 1 1 1 1 1 1 0 1 1 1 @8 CUBIT_ID 2 #
+vertex $97 -1 -1 $-1 $35 $98 #
+straight-curve $-1 -1 -1 $-1 0 5 5 -1 0 0 I I #
+coedge $-1 -1 -1 $-1 $36 $50 $38 $72 reversed $17 $-1 #
+integer_attrib-name_attrib-gen-attrib $-1 -1 $-1 $-1 $37 2 1 1 1 1 1 1 1 1 1 1 1 1 1 0 1 1 1 @8 CUBIT_ID 4 #
+vertex $99 -1 -1 $-1 $62 $100 #
+straight-curve $-1 -1 -1 $-1 0 -5 5 1 0 0 I I #
+edge $101 -1 -1 $-1 $42 -5 $59 5 $38 $102 forward @7 unknown T 5 -5 -5 5 -5 5 #
+edge $103 -1 -1 $-1 $43 -5 $58 5 $39 $104 forward @7 unknown T 5 5 -5 5 5 5 #
+integer_attrib-name_attrib-gen-attrib $-1 -1 $-1 $-1 $42 2 1 1 1 1 1 1 1 1 1 1 1 1 1 0 1 1 1 @8 CUBIT_ID 1 #
+point $-1 -1 -1 $-1 5 -5 5 #
+integer_attrib-name_attrib-gen-attrib $-1 -1 $-1 $-1 $43 2 1 1 1 1 1 1 1 1 1 1 1 1 1 0 1 1 1 @8 CUBIT_ID 2 #
+point $-1 -1 -1 $-1 5 5 5 #
+integer_attrib-name_attrib-gen-attrib $-1 -1 $-1 $-1 $46 2 1 1 1 1 1 1 1 1 1 1 1 1 1 0 1 1 1 @8 CUBIT_ID 6 #
+plane-surface $-1 -1 -1 $-1 5 0 0 -1 0 0 0 -0 1 forward_v I I I I #
+coedge $-1 -1 -1 $-1 $51 $49 $53 $85 reversed $26 $-1 #
+edge $105 -1 -1 $-1 $66 -5 $89 5 $64 $106 forward @7 unknown T -5 5 -5 -5 5 5 #
+integer_attrib-name_attrib-gen-attrib $-1 -1 $-1 $-1 $52 2 1 1 1 1 1 1 1 1 1 1 1 1 1 0 1 1 1 @8 CUBIT_ID 9 #
+vertex $107 -1 -1 $-1 $85 $108 #
+straight-curve $-1 -1 -1 $-1 -5 -5 0 0 0 -1 I I #
+edge $109 -1 -1 $-1 $83 -5 $89 5 $80 $110 forward @7 unknown T -5 -5 -5 -5 5 -5 #
+integer_attrib-name_attrib-gen-attrib $-1 -1 $-1 $-1 $54 2 1 1 1 1 1 1 1 1 1 1 1 1 1 0 1 1 1 @8 CUBIT_ID 6 #
+straight-curve $-1 -1 -1 $-1 0 -5 -5 -1 0 0 I I #
+integer_attrib-name_attrib-gen-attrib $-1 -1 $-1 $-1 $56 2 1 1 1 1 1 1 1 1 1 1 1 1 1 0 1 1 1 @8 CUBIT_ID 8 #
+vertex $111 -1 -1 $-1 $56 $112 #
+straight-curve $-1 -1 -1 $-1 0 5 -5 1 0 0 I I #
+integer_attrib-name_attrib-gen-attrib $-1 -1 $-1 $-1 $58 2 1 1 1 1 1 1 1 1 1 1 1 1 1 0 1 1 1 @8 CUBIT_ID 5 #
+point $-1 -1 -1 $-1 5 5 -5 #
+integer_attrib-name_attrib-gen-attrib $-1 -1 $-1 $-1 $59 2 1 1 1 1 1 1 1 1 1 1 1 1 1 0 1 1 1 @8 CUBIT_ID 6 #
+point $-1 -1 -1 $-1 5 -5 -5 #
+integer_attrib-name_attrib-gen-attrib $-1 -1 $-1 $-1 $62 2 1 1 1 1 1 1 1 1 1 1 1 1 1 0 1 1 1 @8 CUBIT_ID 3 #
+straight-curve $-1 -1 -1 $-1 -5 0 5 0 -1 0 I I #
+integer_attrib-name_attrib-gen-attrib $-1 -1 $-1 $-1 $66 2 1 1 1 1 1 1 1 1 1 1 1 1 1 0 1 1 1 @8 CUBIT_ID 3 #
+point $-1 -1 -1 $-1 -5 5 5 #
+integer_attrib-name_attrib-gen-attrib $-1 -1 $-1 $-1 $70 2 1 1 1 1 1 1 1 1 1 1 1 1 1 0 1 1 1 @8 CUBIT_ID 4 #
+point $-1 -1 -1 $-1 -5 -5 5 #
+integer_attrib-name_attrib-gen-attrib $-1 -1 $-1 $-1 $72 2 1 1 1 1 1 1 1 1 1 1 1 1 1 0 1 1 1 @8 CUBIT_ID 10 #
+straight-curve $-1 -1 -1 $-1 5 -5 0 0 0 -1 I I #
+integer_attrib-name_attrib-gen-attrib $-1 -1 $-1 $-1 $73 2 1 1 1 1 1 1 1 1 1 1 1 1 1 0 1 1 1 @8 CUBIT_ID 12 #
+straight-curve $-1 -1 -1 $-1 5 5 0 0 0 -1 I I #
+integer_attrib-name_attrib-gen-attrib $-1 -1 $-1 $-1 $81 2 1 1 1 1 1 1 1 1 1 1 1 1 1 0 1 1 1 @8 CUBIT_ID 11 #
+straight-curve $-1 -1 -1 $-1 -5 5 0 0 0 -1 I I #
+integer_attrib-name_attrib-gen-attrib $-1 -1 $-1 $-1 $83 2 1 1 1 1 1 1 1 1 1 1 1 1 1 0 1 1 1 @8 CUBIT_ID 7 #
+point $-1 -1 -1 $-1 -5 -5 -5 #
+integer_attrib-name_attrib-gen-attrib $-1 -1 $-1 $-1 $85 2 1 1 1 1 1 1 1 1 1 1 1 1 1 0 1 1 1 @8 CUBIT_ID 7 #
+straight-curve $-1 -1 -1 $-1 -5 0 -5 0 1 0 I I #
+integer_attrib-name_attrib-gen-attrib $-1 -1 $-1 $-1 $89 2 1 1 1 1 1 1 1 1 1 1 1 1 1 0 1 1 1 @8 CUBIT_ID 8 #
+point $-1 -1 -1 $-1 -5 5 -5 #
+End-of-ACIS-data
\ No newline at end of file
diff --git a/MeshFiles/unittest/io/cube.stp b/MeshFiles/unittest/io/cube.stp
new file mode 100644
index 0000000..c17645d
--- /dev/null
+++ b/MeshFiles/unittest/io/cube.stp
@@ -0,0 +1,185 @@
+ISO-10303-21;
+HEADER;
+FILE_DESCRIPTION(('STEP AP214'),'1');
+FILE_NAME('/home/iulian/source/MOABsource/MeshFiles/unittest/io/cube.stp','2014-02-11T15:06:47',(' '),(' '),'Spatial InterOp 3D',' ',' ');
+FILE_SCHEMA(('automotive_design'));
+ENDSEC;
+DATA;
+#1=PRODUCT_DEFINITION_CONTEXT('',#9,'design');
+#2=APPLICATION_PROTOCOL_DEFINITION('INTERNATIONAL STANDARD','automotive_design',1994,#9);
+#3=PRODUCT_CATEGORY_RELATIONSHIP('NONE','NONE',#10,#11);
+#4=SHAPE_DEFINITION_REPRESENTATION(#12,#13);
+#5= (GEOMETRIC_REPRESENTATION_CONTEXT(3)GLOBAL_UNCERTAINTY_ASSIGNED_CONTEXT((#16))GLOBAL_UNIT_ASSIGNED_CONTEXT((#18,#19,#20))REPRESENTATION_CONTEXT('NONE','WORKSPACE'));
+#9=APPLICATION_CONTEXT(' ');
+#10=PRODUCT_CATEGORY('part','NONE');
+#11=PRODUCT_RELATED_PRODUCT_CATEGORY('detail',' ',(#22));
+#12=PRODUCT_DEFINITION_SHAPE('NONE','NONE',#23);
+#13=ADVANCED_BREP_SHAPE_REPRESENTATION('1',(#24,#25),#5);
+#16=UNCERTAINTY_MEASURE_WITH_UNIT(LENGTH_MEASURE(1.0E-06),#18,'','');
+#18= (CONVERSION_BASED_UNIT('MILLIMETRE',#28)LENGTH_UNIT()NAMED_UNIT(#31));
+#19= (NAMED_UNIT(#33)PLANE_ANGLE_UNIT()SI_UNIT($,.RADIAN.));
+#20= (NAMED_UNIT(#33)SOLID_ANGLE_UNIT()SI_UNIT($,.STERADIAN.));
+#22=PRODUCT('1','1','PART-1-DESC',(#39));
+#23=PRODUCT_DEFINITION('NONE','NONE',#40,#1);
+#24=MANIFOLD_SOLID_BREP('1',#41);
+#25=AXIS2_PLACEMENT_3D('',#42,#43,#44);
+#28=LENGTH_MEASURE_WITH_UNIT(LENGTH_MEASURE(1.0),#45);
+#31=DIMENSIONAL_EXPONENTS(1.0,0.0,0.0,0.0,0.0,0.0,0.0);
+#33=DIMENSIONAL_EXPONENTS(0.0,0.0,0.0,0.0,0.0,0.0,0.0);
+#39=PRODUCT_CONTEXT('',#9,'mechanical');
+#40=PRODUCT_DEFINITION_FORMATION_WITH_SPECIFIED_SOURCE(' ','NONE',#22,.NOT_KNOWN.);
+#41=CLOSED_SHELL('',(#46,#47,#48,#49,#50,#51));
+#42=CARTESIAN_POINT('',(0.0,0.0,0.0));
+#43=DIRECTION('',(0.0,0.0,1.0));
+#44=DIRECTION('',(1.0,0.0,0.0));
+#45= (NAMED_UNIT(#31)LENGTH_UNIT()SI_UNIT(.MILLI.,.METRE.));
+#46=ADVANCED_FACE('',(#53),#54,.T.);
+#47=ADVANCED_FACE('',(#55),#56,.F.);
+#48=ADVANCED_FACE('',(#57),#58,.F.);
+#49=ADVANCED_FACE('',(#59),#60,.F.);
+#50=ADVANCED_FACE('',(#61),#62,.F.);
+#51=ADVANCED_FACE('',(#63),#64,.F.);
+#53=FACE_OUTER_BOUND('',#65,.T.);
+#54=PLANE('',#66);
+#55=FACE_OUTER_BOUND('',#67,.T.);
+#56=PLANE('',#68);
+#57=FACE_OUTER_BOUND('',#69,.T.);
+#58=PLANE('',#70);
+#59=FACE_OUTER_BOUND('',#71,.T.);
+#60=PLANE('',#72);
+#61=FACE_OUTER_BOUND('',#73,.T.);
+#62=PLANE('',#74);
+#63=FACE_OUTER_BOUND('',#75,.T.);
+#64=PLANE('',#76);
+#65=EDGE_LOOP('',(#77,#78,#79,#80));
+#66=AXIS2_PLACEMENT_3D('',#81,#82,#83);
+#67=EDGE_LOOP('',(#84,#85,#86,#87));
+#68=AXIS2_PLACEMENT_3D('',#88,#89,#90);
+#69=EDGE_LOOP('',(#91,#92,#93,#94));
+#70=AXIS2_PLACEMENT_3D('',#95,#96,#97);
+#71=EDGE_LOOP('',(#98,#99,#100,#101));
+#72=AXIS2_PLACEMENT_3D('',#102,#103,#104);
+#73=EDGE_LOOP('',(#105,#106,#107,#108));
+#74=AXIS2_PLACEMENT_3D('',#109,#110,#111);
+#75=EDGE_LOOP('',(#112,#113,#114,#115));
+#76=AXIS2_PLACEMENT_3D('',#116,#117,#118);
+#77=ORIENTED_EDGE('',*,*,#119,.T.);
+#78=ORIENTED_EDGE('',*,*,#120,.T.);
+#79=ORIENTED_EDGE('',*,*,#121,.T.);
+#80=ORIENTED_EDGE('',*,*,#122,.T.);
+#81=CARTESIAN_POINT('',(0.0,0.0,5.0));
+#82=DIRECTION('',(0.0,0.0,1.0));
+#83=DIRECTION('',(1.0,0.0,0.0));
+#84=ORIENTED_EDGE('',*,*,#123,.T.);
+#85=ORIENTED_EDGE('',*,*,#124,.T.);
+#86=ORIENTED_EDGE('',*,*,#125,.T.);
+#87=ORIENTED_EDGE('',*,*,#126,.T.);
+#88=CARTESIAN_POINT('',(0.0,0.0,-5.0));
+#89=DIRECTION('',(0.0,0.0,1.0));
+#90=DIRECTION('',(1.0,0.0,0.0));
+#91=ORIENTED_EDGE('',*,*,#127,.T.);
+#92=ORIENTED_EDGE('',*,*,#124,.F.);
+#93=ORIENTED_EDGE('',*,*,#128,.F.);
+#94=ORIENTED_EDGE('',*,*,#122,.F.);
+#95=CARTESIAN_POINT('',(0.0,-5.0,0.0));
+#96=DIRECTION('',(0.0,1.0,-0.0));
+#97=DIRECTION('',(-0.0,0.0,1.0));
+#98=ORIENTED_EDGE('',*,*,#129,.T.);
+#99=ORIENTED_EDGE('',*,*,#125,.F.);
+#100=ORIENTED_EDGE('',*,*,#127,.F.);
+#101=ORIENTED_EDGE('',*,*,#121,.F.);
+#102=CARTESIAN_POINT('',(-5.0,0.0,0.0));
+#103=DIRECTION('',(1.0,0.0,0.0));
+#104=DIRECTION('',(0.0,0.0,-1.0));
+#105=ORIENTED_EDGE('',*,*,#130,.T.);
+#106=ORIENTED_EDGE('',*,*,#126,.F.);
+#107=ORIENTED_EDGE('',*,*,#129,.F.);
+#108=ORIENTED_EDGE('',*,*,#120,.F.);
+#109=CARTESIAN_POINT('',(0.0,5.0,0.0));
+#110=DIRECTION('',(0.0,-1.0,0.0));
+#111=DIRECTION('',(0.0,0.0,-1.0));
+#112=ORIENTED_EDGE('',*,*,#128,.T.);
+#113=ORIENTED_EDGE('',*,*,#123,.F.);
+#114=ORIENTED_EDGE('',*,*,#130,.F.);
+#115=ORIENTED_EDGE('',*,*,#119,.F.);
+#116=CARTESIAN_POINT('',(5.0,0.0,0.0));
+#117=DIRECTION('',(-1.0,0.0,0.0));
+#118=DIRECTION('',(0.0,-0.0,1.0));
+#119=EDGE_CURVE('',#131,#132,#133,.T.);
+#120=EDGE_CURVE('',#132,#134,#135,.T.);
+#121=EDGE_CURVE('',#134,#136,#137,.T.);
+#122=EDGE_CURVE('',#136,#131,#138,.T.);
+#123=EDGE_CURVE('',#139,#140,#141,.T.);
+#124=EDGE_CURVE('',#140,#142,#143,.T.);
+#125=EDGE_CURVE('',#142,#144,#145,.T.);
+#126=EDGE_CURVE('',#144,#139,#146,.T.);
+#127=EDGE_CURVE('',#136,#142,#147,.T.);
+#128=EDGE_CURVE('',#131,#140,#148,.T.);
+#129=EDGE_CURVE('',#134,#144,#149,.T.);
+#130=EDGE_CURVE('',#132,#139,#150,.T.);
+#131=VERTEX_POINT('',#151);
+#132=VERTEX_POINT('',#152);
+#133=LINE('',#153,#154);
+#134=VERTEX_POINT('',#155);
+#135=LINE('',#156,#157);
+#136=VERTEX_POINT('',#158);
+#137=LINE('',#159,#160);
+#138=LINE('',#161,#162);
+#139=VERTEX_POINT('',#163);
+#140=VERTEX_POINT('',#164);
+#141=LINE('',#165,#166);
+#142=VERTEX_POINT('',#167);
+#143=LINE('',#168,#169);
+#144=VERTEX_POINT('',#170);
+#145=LINE('',#171,#172);
+#146=LINE('',#173,#174);
+#147=LINE('',#175,#176);
+#148=LINE('',#177,#178);
+#149=LINE('',#179,#180);
+#150=LINE('',#181,#182);
+#151=CARTESIAN_POINT('',(5.0,-5.0,5.0));
+#152=CARTESIAN_POINT('',(5.0,5.0,5.0));
+#153=CARTESIAN_POINT('',(5.0,0.0,5.0));
+#154=VECTOR('',#183,1.0);
+#155=CARTESIAN_POINT('',(-5.0,5.0,5.0));
+#156=CARTESIAN_POINT('',(0.0,5.0,5.0));
+#157=VECTOR('',#184,1.0);
+#158=CARTESIAN_POINT('',(-5.0,-5.0,5.0));
+#159=CARTESIAN_POINT('',(-5.0,0.0,5.0));
+#160=VECTOR('',#185,1.0);
+#161=CARTESIAN_POINT('',(0.0,-5.0,5.0));
+#162=VECTOR('',#186,1.0);
+#163=CARTESIAN_POINT('',(5.0,5.0,-5.0));
+#164=CARTESIAN_POINT('',(5.0,-5.0,-5.0));
+#165=CARTESIAN_POINT('',(5.0,0.0,-5.0));
+#166=VECTOR('',#187,1.0);
+#167=CARTESIAN_POINT('',(-5.0,-5.0,-5.0));
+#168=CARTESIAN_POINT('',(0.0,-5.0,-5.0));
+#169=VECTOR('',#188,1.0);
+#170=CARTESIAN_POINT('',(-5.0,5.0,-5.0));
+#171=CARTESIAN_POINT('',(-5.0,0.0,-5.0));
+#172=VECTOR('',#189,1.0);
+#173=CARTESIAN_POINT('',(0.0,5.0,-5.0));
+#174=VECTOR('',#190,1.0);
+#175=CARTESIAN_POINT('',(-5.0,-5.0,0.0));
+#176=VECTOR('',#191,1.0);
+#177=CARTESIAN_POINT('',(5.0,-5.0,0.0));
+#178=VECTOR('',#192,1.0);
+#179=CARTESIAN_POINT('',(-5.0,5.0,0.0));
+#180=VECTOR('',#193,1.0);
+#181=CARTESIAN_POINT('',(5.0,5.0,0.0));
+#182=VECTOR('',#194,1.0);
+#183=DIRECTION('',(0.0,1.0,0.0));
+#184=DIRECTION('',(-1.0,0.0,0.0));
+#185=DIRECTION('',(0.0,-1.0,0.0));
+#186=DIRECTION('',(1.0,0.0,0.0));
+#187=DIRECTION('',(0.0,-1.0,0.0));
+#188=DIRECTION('',(-1.0,0.0,0.0));
+#189=DIRECTION('',(0.0,1.0,0.0));
+#190=DIRECTION('',(1.0,0.0,0.0));
+#191=DIRECTION('',(0.0,0.0,-1.0));
+#192=DIRECTION('',(0.0,0.0,-1.0));
+#193=DIRECTION('',(0.0,0.0,-1.0));
+#194=DIRECTION('',(0.0,0.0,-1.0));
+ENDSEC;
+END-ISO-10303-21;
diff --git a/config/compiler.m4 b/config/compiler.m4
index ce83a5f..4787100 100644
--- a/config/compiler.m4
+++ b/config/compiler.m4
@@ -548,6 +548,8 @@ case "$cc_compiler:$host_cpu" in
bgq)
FATHOM_CC_32BIT=-q32
FATHOM_CC_64BIT=-q64
+ FATHOM_CC_SPECIAL=-qarch=qp
+ FATHOM_CXX_SPECIAL="-qarch=qp -qpic=large -qmaxmem=-1"
AR="ar"
NM="nm -B"
;;
diff --git a/examples/makefile b/examples/makefile
index a6e3d1f..6c1ed9a 100644
--- a/examples/makefile
+++ b/examples/makefile
@@ -36,7 +36,7 @@ DirectAccessNoHoles: DirectAccessNoHoles.o ${MOAB_LIBDIR}/libMOAB.la
${MOAB_CXX} -o $@ $< ${MOAB_LIBS_LINK}
DirectAccessNoHolesF90: DirectAccessNoHolesF90.o ${MOAB_LIBDIR}/libMOAB.la
- ${MOAB_CXX} -o $@ $< ${IMESH_LIBS}
+ ${MOAB_FC} -o $@ $< ${IMESH_LIBS}
ReduceExchangeTags : ReduceExchangeTags.o ${MOAB_LIBDIR}/libMOAB.la
${MOAB_CXX} -o $@ $< ${MOAB_LIBS_LINK}
@@ -51,7 +51,7 @@ point_in_elem_search: point_in_elem_search.o ${MOAB_LIBDIR}/libMOAB.la
${MOAB_CXX} -o $@ $< ${MOAB_LIBS_LINK}
PushParMeshIntoMoabF90: PushParMeshIntoMoabF90.o
- ${MOAB_CXX} -o $@ $< ${IMESH_LIBS} -lgfortran -L/usr/lib/openmpi/lib -lmpi_f90 -lmpi_f77 -lmpi -lopen-rte -lopen-pal -ldl
+ ${MOAB_FC} -o $@ $< ${IMESH_LIBS}
DeformMeshRemap: DeformMeshRemap.o ${MOAB_LIBDIR}/libMOAB.la
${MOAB_CXX} -o $@ $< ${MOAB_LIBS_LINK} -lmbcoupler ${MOAB_LIBS_LINK}
diff --git a/examples/point_in_elem_search.cpp b/examples/point_in_elem_search.cpp
index a917fa0..ce22928 100644
--- a/examples/point_in_elem_search.cpp
+++ b/examples/point_in_elem_search.cpp
@@ -27,7 +27,7 @@ int main(int argc, char **argv) {
int num_queries = 1000000;
- if (argc == 1) {
+ if (argc < 2 || argc > 3) {
std::cout << "Usage: " << argv[0] << "<filename> [num_queries]" << std::endl;
return 0;
}
@@ -37,7 +37,7 @@ int main(int argc, char **argv) {
moab::Core mb;
// load the file
- ErrorCode rval = mb.load_file(argv[argc-1]); ERR("Error loading file");
+ ErrorCode rval = mb.load_file(argv[1]); ERR("Error loading file");
// get all 3d elements in the file
Range elems;
@@ -66,7 +66,7 @@ int main(int argc, char **argv) {
for (int i = 0; i < num_queries; i++) {
pos = box.bMin +
CartVect(box_extents[0]*.01*(rand()%100), box_extents[1]*.01*(rand()%100), box_extents[2]*.01*(rand()%100));
- ErrorCode tmp_rval = sl.locate_point(pos.array(), elem, params.array(), 0.0, 0.0, &is_inside);
+ ErrorCode tmp_rval = sl.locate_point(pos.array(), elem, params.array(), &is_inside, 0.0, 0.0);
if (MB_SUCCESS != tmp_rval) rval = tmp_rval;
if (is_inside) num_inside++;
}
diff --git a/src/CMakeLists.txt b/src/CMakeLists.txt
index 280b203..40f7a59 100644
--- a/src/CMakeLists.txt
+++ b/src/CMakeLists.txt
@@ -7,6 +7,7 @@
AxisBox.cpp
BitPage.cpp
BitTag.cpp
+ BoundBox.cpp
BSPTree.cpp
BSPTreePoly.cpp
CN.cpp
@@ -46,6 +47,8 @@
SweptVertexData.cpp
SysUtil.cpp
TagInfo.cpp
+ Tree.cpp
+ TupleList.cpp
Types.cpp
TypeSequenceManager.cpp
UnstructuredElemSeq.cpp
@@ -93,14 +96,27 @@
io/WriteVtk.cpp
ReaderWriterSet.cpp
)
+
+ set ( MOABLD_LIB_SRCS
+ LocalDiscretization/ElemEvaluator.cpp
+ LocalDiscretization/LinearHex.cpp
+ LocalDiscretization/LinearQuad.cpp
+ LocalDiscretization/LinearTet.cpp
+ LocalDiscretization/LinearTri.cpp
+ LocalDiscretization/QuadraticHex.cpp
+# LocalDiscretization/SpectralHex.cpp
+# LocalDiscretization/SpectralQuad.cpp
+ )
+
include_directories(
- ${MOAB_BINARY_DIR}
${MOAB_SOURCE_DIR}/src
- ${MOAB_BINARY_DIR}/src
${MOAB_SOURCE_DIR}/src/io
${MOAB_SOURCE_DIR}/src/parallel
- ${MOAB_BINARY_DIR}/src/parallel
+ ${MOAB_SOURCE_DIR}/src/LocalDiscretization
${MOAB_SOURCE_DIR}/src/moab/point_locater/lotte
+ ${MOAB_BINARY_DIR}
+ ${MOAB_BINARY_DIR}/src
+ ${MOAB_BINARY_DIR}/src/parallel
)
if ( NetCDF_FOUND )
@@ -143,7 +159,7 @@
endif ( HDF5_FOUND )
- SET(MOAB_LIB_SRCS ${MOABCORE_LIB_SRCS} ${MOABPTLOC_LIB_SRCS} ${MOABIO_LIB_SRCS})
+ SET(MOAB_LIB_SRCS ${MOABCORE_LIB_SRCS} ${MOABPTLOC_LIB_SRCS} ${MOABIO_LIB_SRCS} ${MOABLD_LIB_SRCS})
set_source_files_properties( ${MOAB_LIB_SRCS}
COMPILE_FLAGS "-DIS_BUILDING_MB ${MOAB_DEFINES}"
diff --git a/src/Makefile.am b/src/Makefile.am
index 0d474b0..5772d12 100644
--- a/src/Makefile.am
+++ b/src/Makefile.am
@@ -124,6 +124,7 @@ libMOAB_la_SOURCES = \
TagInfo.cpp \
TagInfo.hpp \
Tree.cpp \
+ TupleList.cpp \
Types.cpp \
TypeSequenceManager.cpp \
TypeSequenceManager.hpp \
@@ -195,9 +196,11 @@ nobase_libMOAB_la_include_HEADERS = \
moab/SpectralMeshTool.hpp \
moab/Tree.hpp \
moab/TreeStats.hpp \
+ moab/TupleList.hpp \
moab/Types.hpp \
moab/UnknownInterface.hpp \
moab/Util.hpp \
+ moab/Version.h \
moab/WriteUtilIface.hpp \
moab/WriterIface.hpp \
MBEntityType.h \
diff --git a/src/SpatialLocator.cpp b/src/SpatialLocator.cpp
index ae124dc..694122c 100644
--- a/src/SpatialLocator.cpp
+++ b/src/SpatialLocator.cpp
@@ -319,7 +319,7 @@ namespace moab
locTable.enableWriteAccess();
// pass storage directly into locate_points, since we know those arrays are contiguous
- ErrorCode rval = locate_points(pos, num_points, locTable.vul_wr, locTable.vr_wr, NULL, rel_iter_tol, abs_iter_tol,
+ ErrorCode rval = locate_points(pos, num_points, (EntityHandle*)locTable.vul_wr, locTable.vr_wr, NULL, rel_iter_tol, abs_iter_tol,
inside_tol);
std::fill(locTable.vi_wr, locTable.vi_wr+num_points, 0);
locTable.set_n(num_points);
diff --git a/src/TupleList.cpp b/src/TupleList.cpp
new file mode 100644
index 0000000..bf08946
--- /dev/null
+++ b/src/TupleList.cpp
@@ -0,0 +1,860 @@
+#include <string.h>
+#include <limits.h>
+#include <stdio.h>
+#include <stdlib.h>
+#include <stdarg.h>
+#include <iostream>
+
+#include "moab/TupleList.hpp"
+
+namespace moab {
+
+void fail(const char *fmt, ...)
+{
+ va_list ap;
+ va_start(ap, fmt);
+ vfprintf(stderr, fmt, ap);
+ va_end(ap);
+ exit(1);
+}
+
+TupleList::buffer::buffer(size_t sz)
+{
+ ptr = NULL;
+ buffSize = 0;
+ this->buffer_init_(sz, __FILE__);
+}
+
+TupleList::buffer::buffer()
+{
+ buffSize = 0;
+ ptr = NULL;
+}
+
+void TupleList::buffer::buffer_init_(size_t sizeIn, const char *file)
+{
+ this->buffSize = sizeIn;
+ void *res = malloc(this->buffSize);
+ if (!res && buffSize > 0)
+ fail("%s: allocation of %d bytes failed\n", file, (int) buffSize);
+ ptr = (char*) res;
+}
+
+void TupleList::buffer::buffer_reserve_(size_t min, const char *file)
+{
+ if (this->buffSize < min)
+ {
+ size_t newSize = this->buffSize;
+ newSize += newSize / 2 + 1;
+ if (newSize < min)
+ newSize = min;
+ void *res = realloc(ptr, newSize);
+ if (!res && newSize > 0)
+ fail("%s: allocation of %d bytes failed\n", file, (int) this->buffSize);
+ ptr = (char*) res;
+ this->buffSize = newSize;
+ }
+}
+
+void TupleList::buffer::reset()
+{
+ free(ptr);
+ ptr = NULL;
+ buffSize = 0;
+}
+
+TupleList::TupleList(uint p_mi, uint p_ml, uint p_mul, uint p_mr, uint p_max)
+{
+ vi = NULL;
+ vl = NULL;
+ vul = NULL;
+ vr = NULL;
+ initialize(p_mi, p_ml, p_mul, p_mr, p_max);
+}
+
+TupleList::TupleList()
+{
+ vi = NULL;
+ vl = NULL;
+ vul = NULL;
+ vr = NULL;
+ disableWriteAccess();
+}
+
+// Allocates space for the tuple list in memory according to parameters
+void TupleList::initialize(uint p_mi, uint p_ml, uint p_mul, uint p_mr,
+ uint p_max)
+{
+ this->n = 0;
+ this->max = p_max;
+ this->mi = p_mi;
+ this->ml = p_ml;
+ this->mul = p_mul;
+ this->mr = p_mr;
+ size_t sz;
+
+ if (max * mi > 0)
+ {
+ sz = max * mi * sizeof(sint);
+ void *resi = malloc(sz);
+ if (!resi && max * mi > 0)
+ fail("%s: allocation of %d bytes failed\n", __FILE__, (int) sz);
+ vi = (sint*) resi;
+ }
+ else
+ vi = NULL;
+ if (max * ml > 0)
+ {
+ sz = max * ml * sizeof(slong);
+ void *resl = malloc(sz);
+ if (!resl && max * ml > 0)
+ fail("%s: allocation of %d bytes failed\n", __FILE__, (int) sz);
+ vl = (slong*) resl;
+ }
+ else
+ vl = NULL;
+ if (max * mul > 0)
+ {
+ sz = max * mul * sizeof(ulong);
+ void *resu = malloc(sz);
+ if (!resu && max * mul > 0)
+ fail("%s: allocation of %d bytes failed\n", __FILE__, (int) sz);
+ vul = (ulong*) resu;
+ }
+ else
+ vul = NULL;
+ if (max * mr > 0)
+ {
+ sz = max * mr * sizeof(realType);
+ void *resr = malloc(sz);
+ if (!resr && max * ml > 0)
+ fail("%s: allocation of %d bytes failed\n", __FILE__, (int) sz);
+ vr = (realType*) resr;
+ }
+ else
+ vr = NULL;
+
+ //Begin with write access disabled
+ this->disableWriteAccess();
+
+ //Set read variables
+ vi_rd = vi;
+ vl_rd = vl;
+ vul_rd = vul;
+ vr_rd = vr;
+}
+
+// Resizes a tuplelist to the given uint max
+ErrorCode TupleList::resize(uint maxIn)
+{
+ this->max = maxIn;
+ size_t sz;
+
+ if (vi || (max * mi > 0))
+ {
+ sz = max * mi * sizeof(sint);
+ void *resi = realloc(vi, sz);
+ if (!resi && max * mi > 0)
+ {
+ fail("%s: allocation of %d bytes failed\n", __FILE__, (int) sz);
+ return moab::MB_MEMORY_ALLOCATION_FAILED;
+ }
+ vi = (sint*) resi;
+ }
+ if (vl || (max * ml > 0))
+ {
+ sz = max * ml * sizeof(slong);
+ void *resl = realloc(vl, sz);
+ if (!resl && max * ml > 0)
+ {
+ fail("%s: allocation of %d bytes failed\n", __FILE__, (int) sz);
+ return moab::MB_MEMORY_ALLOCATION_FAILED;
+ }
+ vl = (slong*) resl;
+ }
+ if (vul || (max * mul > 0))
+ {
+ sz = max * mul * sizeof(ulong);
+ void *resu = realloc(vul, sz);
+ if (!resu && max * mul > 0)
+ {
+ fail("%s: allocation of %d bytes failed\n", __FILE__, (int) sz);
+ return moab::MB_MEMORY_ALLOCATION_FAILED;
+ }
+ vul = (ulong*) resu;
+ }
+ if (vr || (max * mr > 0))
+ {
+ sz = max * mr * sizeof(realType);
+ void *resr = realloc(vr, sz);
+ if (!resr && max * mr > 0)
+ {
+ fail("%s: allocation of %d bytes failed\n", __FILE__, (int) sz);
+ return moab::MB_MEMORY_ALLOCATION_FAILED;
+ }
+ vr = (realType*) resr;
+ }
+
+ //Set read variables
+ vi_rd = vi;
+ vl_rd = vl;
+ vul_rd = vul;
+ vr_rd = vr;
+
+ //Set the write variables if necessary
+ if (writeEnabled)
+ {
+ vi_wr = vi;
+ vl_wr = vl;
+ vul_wr = vul;
+ vr_wr = vr;
+ }
+ return moab::MB_SUCCESS;
+}
+
+// Frees the memory used by the tuplelist
+void TupleList::reset()
+{
+ //free up the pointers
+ free(vi);
+ free(vl);
+ free(vul);
+ free(vr);
+ //Set them all to null
+ vr = NULL;
+ vi = NULL;
+ vul = NULL;
+ vl = NULL;
+ //Set the read and write pointers to null
+ disableWriteAccess();
+ vi_rd = NULL;
+ vl_rd = NULL;
+ vul_rd = NULL;
+ vr_rd = NULL;
+}
+
+// Increments n; if n>max, increase the size of the tuplelist
+void TupleList::reserve()
+{
+ n++;
+ while (n > max)
+ resize((max ? max + max / 2 + 1 : 2));
+ last_sorted = -1;
+}
+
+// Given the value and the position in the field, finds the index of the tuple
+// to which the value belongs
+int TupleList::find(unsigned int key_num, sint value)
+{
+ ulong uvalue = (ulong) value;
+ if (!(key_num > mi))
+ {
+ // Binary search: only if the tuple_list is sorted
+ if (last_sorted == (int) key_num)
+ {
+ int lb = 0, ub = n, index; // lb=lower bound, ub=upper bound, index=mid
+ for (; lb <= ub;)
+ {
+ index = (lb + ub) / 2;
+ if (vi[index * mi + key_num] == (long) uvalue)
+ return index;
+ else if (vi[index * mi + key_num] > (long) uvalue)
+ ub = index - 1;
+ else if (vi[index * mi + key_num] < (long) uvalue)
+ lb = index + 1;
+ }
+ }
+ else
+ {
+ // Sequential search: if tuple_list is not sorted
+ for (long index = 0; index < n; index++)
+ {
+ if (vi[index * mi + key_num] == (long) uvalue)
+ return index;
+ }
+ }
+ }
+ return -1; // If the value wasn't present or an invalid key was given
+}
+
+int TupleList::find(unsigned int key_num, slong value)
+{
+ ulong uvalue = (ulong) value;
+ if (!(key_num > ml))
+ {
+ if (last_sorted - mi == key_num)
+ {
+ int lb = 0, ub = n, index; // lb=lower bound, ub=upper bound, index=mid
+ for (; lb <= ub;)
+ {
+ index = (lb + ub) / 2;
+ if (vl[index * ml + key_num] == (long) uvalue)
+ return index;
+ else if (vl[index * ml + key_num] > (long) uvalue)
+ ub = index - 1;
+ else if (vl[index * ml + key_num] < (long) uvalue)
+ lb = index + 1;
+ }
+ }
+ else
+ {
+ // Sequential search: if tuple_list is not sorted
+ for (uint index = 0; index < n; index++)
+ {
+ if (vl[index * ml + key_num] == (long) uvalue)
+ return index;
+ }
+ }
+ }
+ return -1; // If the value wasn't present or an invalid key was given
+}
+
+int TupleList::find(unsigned int key_num, ulong value)
+{
+ if (!(key_num > mul))
+ {
+ if (last_sorted - mi - ml == key_num)
+ {
+ int lb = 0, ub = n - 1, index; // lb=lower bound, ub=upper bound, index=mid
+ for (; lb <= ub;)
+ {
+ index = (lb + ub) / 2;
+ if (vul[index * mul + key_num] == value)
+ return index;
+ else if (vul[index * mul + key_num] > value)
+ ub = index - 1;
+ else if (vul[index * mul + key_num] < value)
+ lb = index + 1;
+ }
+ }
+ else
+ {
+ // Sequential search: if tuple_list is not sorted
+ for (uint index = 0; index < n; index++)
+ {
+ if (vul[index * mul + key_num] == value)
+ return index;
+ }
+ }
+ }
+ return -1; // If the value wasn't present or an invalid key was given
+}
+
+int TupleList::find(unsigned int key_num, realType value)
+{
+ if (!key_num > mr)
+ {
+ // Sequential search: TupleList cannot be sorted by reals
+ for (uint index = 0; index < n; index++)
+ {
+ if (vr[index * mr + key_num] == value)
+ return index;
+ }
+ }
+ return -1; // If the value wasn't present or an invalid key was given
+}
+
+sint TupleList::get_sint(unsigned int index, unsigned int m)
+{
+ if (mi > m && n > index)
+ return vi[index * mi + m];
+ return 0;
+}
+
+slong TupleList::get_int(unsigned int index, unsigned int m)
+{
+ if (ml > m && n > index)
+ return vl[index * ml + m];
+ return 0;
+}
+
+ulong TupleList::get_ulong(unsigned int index, unsigned int m)
+{
+ if (mul > m && n > index)
+ return vul[index * mul + m];
+ return 0;
+}
+
+realType TupleList::get_double(unsigned int index, unsigned int m)
+{
+ if (mr > m && n > index)
+ return vr[index * mr + m];
+ return 0;
+}
+
+ErrorCode TupleList::get(unsigned int index, const sint *&sp, const slong *&ip,
+ const ulong *&lp, const realType *&dp)
+{
+ if (index <= n)
+ {
+ if (mi)
+ *&sp = &vi[index * mi];
+ else
+ *&sp = NULL;
+ if (ml)
+ *&ip = &vl[index * ml];
+ else
+ *&ip = NULL;
+ if (mul)
+ *&lp = &vul[index * mul];
+ else
+ *&lp = NULL;
+ if (mr)
+ *&dp = &vr[index * mr];
+ else
+ *&dp = NULL;
+
+ return MB_SUCCESS;
+ }
+ return MB_FAILURE;
+}
+
+unsigned int TupleList::push_back(sint *sp, slong *ip, ulong *lp, realType *dp)
+{
+ reserve();
+ if (mi)
+ memcpy(&vi[mi * (n - 1)], sp, mi * sizeof(sint));
+ if (ml)
+ memcpy(&vl[ml * (n - 1)], ip, ml * sizeof(long));
+ if (mul)
+ memcpy(&vul[mul * (n - 1)], lp, mul * sizeof(ulong));
+ if (mr)
+ memcpy(&vr[mr * (n - 1)], dp, mr * sizeof(realType));
+
+ last_sorted = -1;
+ return n - 1;
+}
+
+void TupleList::enableWriteAccess()
+{
+ writeEnabled = true;
+ last_sorted = -1;
+ vi_wr = vi;
+ vl_wr = vl;
+ vul_wr = vul;
+ vr_wr = vr;
+}
+
+void TupleList::disableWriteAccess()
+{
+ writeEnabled = false;
+ vi_wr = NULL;
+ vl_wr = NULL;
+ vul_wr = NULL;
+ vr_wr = NULL;
+}
+
+void TupleList::getTupleSize(uint &mi_out, uint &ml_out, uint &mul_out,
+ uint &mr_out) const
+{
+ mi_out = mi;
+ ml_out = ml;
+ mul_out = mul;
+ mr_out = mr;
+}
+
+uint TupleList::inc_n()
+{
+ //Check for direct write access
+ if (!writeEnabled)
+ {
+ enableWriteAccess();
+ }
+ n++;
+ return n;
+}
+
+void TupleList::set_n(uint n_in)
+{
+ //Check for direct write access;
+ if (!writeEnabled)
+ {
+ enableWriteAccess();
+ }
+ n = n_in;
+}
+
+void TupleList::print(const char *name) const
+{
+ std::cout << "Printing Tuple " << name << "===================" << std::endl;
+ unsigned long i = 0, l = 0, ul = 0, r = 0;
+ for (uint k = 0; k < n; k++)
+ {
+ for (uint j = 0; j < mi; j++)
+ {
+ std::cout << vi[i++] << " | ";
+ }
+ for (uint j = 0; j < ml; j++)
+ {
+ std::cout << vl[l++] << " | ";
+ }
+ for (uint j = 0; j < mul; j++)
+ {
+ std::cout << vul[ul++] << " | ";
+ }
+ for (uint j = 0; j < mr; j++)
+ {
+ std::cout << vr[r++] << " | ";
+ }
+ std::cout << std::endl;
+ }
+ std::cout << "=======================================" << std::endl
+ << std::endl;
+}
+
+void TupleList::permute(uint *perm, void *work)
+{
+ const unsigned int_size = mi * sizeof(sint), long_size = ml * sizeof(slong),
+ ulong_size = mul * sizeof(ulong), real_size = mr * sizeof(realType);
+ if (mi)
+ {
+ uint *p = perm, *pe = p + n;
+ char *sorted = (char *) work;
+ while (p != pe)
+ memcpy((void *) sorted, &vi[mi * (*p++)], int_size), sorted += int_size;
+ memcpy(vi, work, int_size * n);
+ }
+ if (ml)
+ {
+ uint *p = perm, *pe = p + n;
+ char *sorted = (char *) work;
+ while (p != pe)
+ memcpy((void *) sorted, &vl[ml * (*p++)], long_size), sorted += long_size;
+ memcpy(vl, work, long_size * n);
+ }
+ if (mul)
+ {
+ uint *p = perm, *pe = p + n;
+ char *sorted = (char *) work;
+ while (p != pe)
+ memcpy((void *) sorted, &vul[mul * (*p++)], ulong_size), sorted +=
+ ulong_size;
+ memcpy(vul, work, ulong_size * n);
+ }
+ if (mr)
+ {
+ uint *p = perm, *pe = p + n;
+ char *sorted = (char *) work;
+ while (p != pe)
+ memcpy((void *) sorted, &vr[mr * (*p++)], real_size), sorted += real_size;
+ memcpy(vr, work, real_size * n);
+ }
+}
+
+# define umax_2(a, b) (((a)>(b)) ? (a):(b))
+
+ErrorCode TupleList::sort(uint key, TupleList::buffer *buf)
+{
+ const unsigned int_size = mi * sizeof(sint);
+ const unsigned long_size = ml * sizeof(slong);
+ const unsigned ulong_size = mul * sizeof(ulong);
+ const unsigned real_size = mr * sizeof(realType);
+ const unsigned width = umax_2(umax_2(int_size,long_size),
+ umax_2(ulong_size,real_size));
+ const unsigned data_size =
+ key >= mi ? sizeof(SortData<long> ) : sizeof(SortData<uint> );
+
+ uint work_min = n * umax_2(2*data_size,sizeof(sint)+width);
+ uint *work;
+ buf->buffer_reserve(work_min);
+ work = (uint *) buf->ptr;
+ if (key < mi)
+ index_sort((uint *) &vi[key], n, mi, work, (SortData<uint>*) work);
+ else if (key < mi + ml)
+ index_sort((long*) &vl[key - mi], n, ml, work, (SortData<long>*) work);
+ else if (key < mi + ml + mul)
+ index_sort((ulong*) &vul[key - mi - ml], n, mul, work,
+ (SortData<ulong>*) work);
+ else
+ return MB_NOT_IMPLEMENTED;
+
+ permute(work, work + n);
+
+ if (!writeEnabled)
+ last_sorted = key;
+ return MB_SUCCESS;
+}
+
+#undef umax_2
+
+#define DIGIT_BITS 8
+#define DIGIT_VALUES (1<<DIGIT_BITS)
+#define DIGIT_MASK ((Value)(DIGIT_VALUES-1))
+#define CEILDIV(a,b) (((a)+(b)-1)/(b))
+#define DIGITS CEILDIV(CHAR_BIT*sizeof(Value),DIGIT_BITS)
+#define VALUE_BITS (DIGIT_BITS*DIGITS)
+#define COUNT_SIZE (DIGITS*DIGIT_VALUES)
+
+/* used to unroll a tiny loop: */
+#define COUNT_DIGIT_01(n,i) \
+ if(n>i) count[i][val&DIGIT_MASK]++, val>>=DIGIT_BITS
+#define COUNT_DIGIT_02(n,i) COUNT_DIGIT_01(n,i); COUNT_DIGIT_01(n,i+ 1)
+#define COUNT_DIGIT_04(n,i) COUNT_DIGIT_02(n,i); COUNT_DIGIT_02(n,i+ 2)
+#define COUNT_DIGIT_08(n,i) COUNT_DIGIT_04(n,i); COUNT_DIGIT_04(n,i+ 4)
+#define COUNT_DIGIT_16(n,i) COUNT_DIGIT_08(n,i); COUNT_DIGIT_08(n,i+ 8)
+#define COUNT_DIGIT_32(n,i) COUNT_DIGIT_16(n,i); COUNT_DIGIT_16(n,i+16)
+#define COUNT_DIGIT_64(n,i) COUNT_DIGIT_32(n,i); COUNT_DIGIT_32(n,i+32)
+
+template<class Value>
+Value TupleList::radix_count(const Value *A, const Value *end, Index stride,
+ Index count[DIGITS][DIGIT_VALUES])
+{
+ Value bitorkey = 0;
+ memset(count, 0, COUNT_SIZE * sizeof(Index));
+ do
+ {
+ Value val = *A;
+ bitorkey |= val;
+ COUNT_DIGIT_64(DIGITS, 0);
+ // above macro expands to:
+ //if(DIGITS> 0) count[ 0][val&DIGIT_MASK]++, val>>=DIGIT_BITS;
+ //if(DIGITS> 1) count[ 1][val&DIGIT_MASK]++, val>>=DIGIT_BITS;
+ // ...
+ //if(DIGITS>63) count[63][val&DIGIT_MASK]++, val>>=DIGIT_BITS;
+
+ } while (A += stride, A != end);
+ return bitorkey;
+}
+
+#undef COUNT_DIGIT_01
+#undef COUNT_DIGIT_02
+#undef COUNT_DIGIT_04
+#undef COUNT_DIGIT_08
+#undef COUNT_DIGIT_16
+#undef COUNT_DIGIT_32
+#undef COUNT_DIGIT_64
+
+void TupleList::radix_offsets(Index *c)
+{
+ Index sum = 0, t, *ce = c + DIGIT_VALUES;
+ do
+ t = *c, *c++ = sum, sum += t;
+ while (c != ce);
+}
+
+template<class Value>
+unsigned TupleList::radix_zeros(Value bitorkey,
+ Index count[DIGITS][DIGIT_VALUES], unsigned *shift, Index **offsets)
+{
+ unsigned digits = 0, sh = 0;
+ Index *c = &count[0][0];
+ do
+ {
+ if (bitorkey & DIGIT_MASK)
+ *shift++ = sh, *offsets++ = c, ++digits, radix_offsets(c);
+ } while (bitorkey >>= DIGIT_BITS,sh += DIGIT_BITS,c += DIGIT_VALUES,sh
+ != VALUE_BITS);
+ return digits;
+}
+
+template<class Value>
+void TupleList::radix_index_pass_b(const Value *A, Index n, Index stride,
+ unsigned sh, Index *off, SortData<Value> *out)
+{
+ Index i = 0;
+ do
+ {
+ Value v = *A;
+ SortData<Value> *d = &out[off[(v >> sh) & DIGIT_MASK]++];
+ d->v = v, d->i = i++;
+ } while (A += stride, i != n);
+}
+
+template<class Value>
+void TupleList::radix_index_pass_m(const SortData<Value> *src,
+ const SortData<Value> *end, unsigned sh, Index *off, SortData<Value> *out)
+{
+ do
+ {
+ SortData<Value> *d = &out[off[(src->v >> sh) & DIGIT_MASK]++];
+ d->v = src->v, d->i = src->i;
+ } while (++src != end);
+}
+
+template<class Value>
+void TupleList::radix_index_pass_e(const SortData<Value> *src,
+ const SortData<Value> *end, unsigned sh, Index *off, Index *out)
+{
+ do
+ out[off[(src->v >> sh) & DIGIT_MASK]++] = src->i;
+ while (++src != end);
+}
+
+template<class Value>
+void TupleList::radix_index_pass_be(const Value *A, Index n, Index stride,
+ unsigned sh, Index *off, Index *out)
+{
+ Index i = 0;
+ do
+ out[off[(*A >> sh) & DIGIT_MASK]++] = i++;
+ while (A += stride, i != n);
+}
+
+template<class Value>
+void TupleList::radix_index_sort(const Value *A, Index n, Index stride,
+ Index *idx, SortData<Value> *work)
+{
+ Index count[DIGITS][DIGIT_VALUES];
+ Value bitorkey = radix_count(A, A + n * stride, stride, count);
+ unsigned shift[DIGITS];
+ Index *offsets[DIGITS];
+ unsigned digits = radix_zeros(bitorkey, count, shift, offsets);
+ if (digits == 0)
+ {
+ Index i = 0;
+ do
+ *idx++ = i++;
+ while (i != n);
+ }
+ else if (digits == 1)
+ {
+ radix_index_pass_be(A, n, stride, shift[0], offsets[0], idx);
+ }
+ else
+ {
+ SortData<Value> *src, *dst;
+ unsigned d;
+ if ((digits & 1) == 0)
+ dst = work, src = dst + n;
+ else
+ src = work, dst = src + n;
+ radix_index_pass_b(A, n, stride, shift[0], offsets[0], src);
+ for (d = 1; d != digits - 1; ++d)
+ {
+ SortData<Value> *t;
+ radix_index_pass_m(src, src + n, shift[d], offsets[d], dst);
+ t = src, src = dst, dst = t;
+ }
+ radix_index_pass_e(src, src + n, shift[d], offsets[d], idx);
+ }
+}
+
+template<class Value>
+void TupleList::merge_index_sort(const Value *A, const Index An, Index stride,
+ Index *idx, SortData<Value> *work)
+{
+ SortData<Value> * const buf[2] = { work + An, work };
+ Index n = An, base = -n, odd = 0, c = 0, b = 1;
+ Index i = 0;
+ for (;;)
+ {
+ SortData<Value> *p;
+ if ((c & 1) == 0)
+ {
+ base += n, n += (odd & 1), c |= 1, b ^= 1;
+ while (n > 3)
+ odd <<= 1, odd |= (n & 1), n >>= 1, c <<= 1, b ^= 1;
+ }
+ else
+ base -= n - (odd & 1), n <<= 1, n -= (odd & 1), odd >>= 1, c >>= 1;
+ if (c == 0)
+ break;
+ p = buf[b] + base;
+ if (n == 2)
+ {
+ Value v[2];
+ v[0] = *A, A += stride, v[1] = *A, A += stride;
+ if (v[1] < v[0])
+ p[0].v = v[1], p[0].i = i + 1, p[1].v = v[0], p[1].i = i;
+ else
+ p[0].v = v[0], p[0].i = i, p[1].v = v[1], p[1].i = i + 1;
+ i += 2;
+ }
+ else if (n == 3)
+ {
+ Value v[3];
+ v[0] = *A, A += stride, v[1] = *A, A += stride, v[2] = *A, A += stride;
+ if (v[1] < v[0])
+ {
+ if (v[2] < v[1])
+ p[0].v = v[2], p[1].v = v[1], p[2].v = v[0], p[0].i = i + 2, p[1].i =
+ i + 1, p[2].i = i;
+ else
+ {
+ if (v[2] < v[0])
+ p[0].v = v[1], p[1].v = v[2], p[2].v = v[0], p[0].i = i + 1, p[1].i =
+ i + 2, p[2].i = i;
+ else
+ p[0].v = v[1], p[1].v = v[0], p[2].v = v[2], p[0].i = i + 1, p[1].i =
+ i, p[2].i = i + 2;
+ }
+ }
+ else
+ {
+ if (v[2] < v[0])
+ p[0].v = v[2], p[1].v = v[0], p[2].v = v[1], p[0].i = i + 2, p[1].i =
+ i, p[2].i = i + 1;
+ else
+ {
+ if (v[2] < v[1])
+ p[0].v = v[0], p[1].v = v[2], p[2].v = v[1], p[0].i = i, p[1].i = i
+ + 2, p[2].i = i + 1;
+ else
+ p[0].v = v[0], p[1].v = v[1], p[2].v = v[2], p[0].i = i, p[1].i = i
+ + 1, p[2].i = i + 2;
+ }
+ }
+ i += 3;
+ }
+ else
+ {
+ const Index na = n >> 1, nb = (n + 1) >> 1;
+ const SortData<Value> *ap = buf[b ^ 1] + base, *ae = ap + na;
+ SortData<Value> *bp = p + na, *be = bp + nb;
+ for (;;)
+ {
+ if (bp->v < ap->v)
+ {
+ *p++ = *bp++;
+ if (bp != be)
+ continue;
+ do
+ *p++ = *ap++;
+ while (ap != ae);
+ break;
+ }
+ else
+ {
+ *p++ = *ap++;
+ if (ap == ae)
+ break;
+ }
+ }
+ }
+ }
+ {
+ const SortData<Value> *p = buf[0], *pe = p + An;
+ do
+ *idx++ = (p++)->i;
+ while (p != pe);
+ }
+}
+
+template<class Value>
+void TupleList::index_sort(const Value *A, Index n, Index stride, Index *idx,
+ SortData<Value> *work)
+{
+ if (n < DIGIT_VALUES)
+ {
+ if (n == 0)
+ return;
+ if (n == 1)
+ *idx = 0;
+ else
+ merge_index_sort(A, n, stride, idx, work);
+ }
+ else
+ radix_index_sort(A, n, stride, idx, work);
+}
+
+#undef DIGIT_BITS
+#undef DIGIT_VALUES
+#undef DIGIT_MASK
+#undef CEILDIV
+#undef DIGITS
+#undef VALUE_BITS
+#undef COUNT_SIZE
+#undef sort_data_long
+
+} //namespace
+
diff --git a/src/io/NCHelper.cpp b/src/io/NCHelper.cpp
index c915382..9ef1332 100644
--- a/src/io/NCHelper.cpp
+++ b/src/io/NCHelper.cpp
@@ -48,6 +48,9 @@ NCHelper* NCHelper::get_nc_helper(ReadNC* readNC, int fileId, const FileOptions&
else {
if (NCHelperMPAS::can_read_file(readNC))
return new (std::nothrow) NCHelperMPAS(readNC, fileId, opts, fileSet);
+ // For a HOMME connectivity file, there might be no CF convention
+ else if (NCHelperHOMME::can_read_file(readNC, fileId))
+ return new (std::nothrow) NCHelperHOMME(readNC, fileId, opts, fileSet);
}
// Unknown NetCDF grid (will fill this in later for POP, CICE and CLM)
diff --git a/src/io/NCHelperHOMME.cpp b/src/io/NCHelperHOMME.cpp
index efddb0c..53fa9ad 100644
--- a/src/io/NCHelperHOMME.cpp
+++ b/src/io/NCHelperHOMME.cpp
@@ -15,7 +15,7 @@ namespace moab {
NCHelperHOMME::NCHelperHOMME(ReadNC* readNC, int fileId, const FileOptions& opts, EntityHandle fileSet)
: UcdNCHelper(readNC, fileId, opts, fileSet),
-_spectralOrder(-1), connectId(-1)
+_spectralOrder(-1), connectId(-1), isConnFile(false)
{
// Calculate spectral order
std::map<std::string, ReadNC::AttData>::iterator attIt = readNC->globalAtts.find("np");
@@ -26,6 +26,11 @@ _spectralOrder(-1), connectId(-1)
else
_spectralOrder--; // Spectral order is one less than np
}
+ else {
+ // As can_read_file() returns true and there is no global attribute "np", it should be a connectivity file
+ isConnFile = true;
+ _spectralOrder = 3; // Assume np is 4
+ }
}
bool NCHelperHOMME::can_read_file(ReadNC* readNC, int fileId)
@@ -52,6 +57,15 @@ bool NCHelperHOMME::can_read_file(ReadNC* readNC, int fileId)
return true;
}
+ else {
+ // If dimension names "ncol" AND "ncorners" AND "ncells" exist, then it should be the HOMME connectivity file
+ // In this case, the mesh can still be created
+ std::vector<std::string>& dimNames = readNC->dimNames;
+ if ((std::find(dimNames.begin(), dimNames.end(), std::string("ncol")) != dimNames.end()) &&
+ (std::find(dimNames.begin(), dimNames.end(), std::string("ncorners")) != dimNames.end()) &&
+ (std::find(dimNames.begin(), dimNames.end(), std::string("ncells")) != dimNames.end()))
+ return true;
+ }
return false;
}
@@ -67,15 +81,20 @@ ErrorCode NCHelperHOMME::init_mesh_vals()
std::vector<std::string>::iterator vit;
// Look for time dimension
- if ((vit = std::find(dimNames.begin(), dimNames.end(), "time")) != dimNames.end())
- idx = vit - dimNames.begin();
- else if ((vit = std::find(dimNames.begin(), dimNames.end(), "t")) != dimNames.end())
- idx = vit - dimNames.begin();
+ if (isConnFile) {
+ // Connectivity file might not have time dimension
+ }
else {
- ERRORR(MB_FAILURE, "Couldn't find 'time' or 't' dimension.");
+ if ((vit = std::find(dimNames.begin(), dimNames.end(), "time")) != dimNames.end())
+ idx = vit - dimNames.begin();
+ else if ((vit = std::find(dimNames.begin(), dimNames.end(), "t")) != dimNames.end())
+ idx = vit - dimNames.begin();
+ else {
+ ERRORR(MB_FAILURE, "Couldn't find 'time' or 't' dimension.");
+ }
+ tDim = idx;
+ nTimeSteps = dimLens[idx];
}
- tDim = idx;
- nTimeSteps = dimLens[idx];
// Get number of vertices (labeled as number of columns)
if ((vit = std::find(dimNames.begin(), dimNames.end(), "ncol")) != dimNames.end())
@@ -90,15 +109,20 @@ ErrorCode NCHelperHOMME::init_mesh_vals()
nCells = nVertices - 2;
// Get number of levels
- if ((vit = std::find(dimNames.begin(), dimNames.end(), "lev")) != dimNames.end())
- idx = vit - dimNames.begin();
- else if ((vit = std::find(dimNames.begin(), dimNames.end(), "ilev")) != dimNames.end())
- idx = vit - dimNames.begin();
+ if (isConnFile) {
+ // Connectivity file might not have level dimension
+ }
else {
- ERRORR(MB_FAILURE, "Couldn't find 'lev' or 'ilev' dimension.");
+ if ((vit = std::find(dimNames.begin(), dimNames.end(), "lev")) != dimNames.end())
+ idx = vit - dimNames.begin();
+ else if ((vit = std::find(dimNames.begin(), dimNames.end(), "ilev")) != dimNames.end())
+ idx = vit - dimNames.begin();
+ else {
+ ERRORR(MB_FAILURE, "Couldn't find 'lev' or 'ilev' dimension.");
+ }
+ levDim = idx;
+ nLevels = dimLens[idx];
}
- levDim = idx;
- nLevels = dimLens[idx];
// Store lon values in xVertVals
std::map<std::string, ReadNC::VarData>::iterator vmit;
@@ -120,35 +144,45 @@ ErrorCode NCHelperHOMME::init_mesh_vals()
}
// Store lev values in levVals
- if ((vmit = varInfo.find("lev")) != varInfo.end() && (*vmit).second.varDims.size() == 1) {
- rval = read_coordinate("lev", 0, nLevels - 1, levVals);
- ERRORR(rval, "Trouble reading 'lev' variable.");
-
- // Decide whether down is positive
- char posval[10];
- int success = NCFUNC(get_att_text)(_fileId, (*vmit).second.varId, "positive", posval);
- if (0 == success && !strcmp(posval, "down")) {
- for (std::vector<double>::iterator dvit = levVals.begin(); dvit != levVals.end(); ++dvit)
- (*dvit) *= -1.0;
- }
+ if (isConnFile) {
+ // Connectivity file might not have level variable
}
else {
- ERRORR(MB_FAILURE, "Couldn't find 'lev' variable.");
+ if ((vmit = varInfo.find("lev")) != varInfo.end() && (*vmit).second.varDims.size() == 1) {
+ rval = read_coordinate("lev", 0, nLevels - 1, levVals);
+ ERRORR(rval, "Trouble reading 'lev' variable.");
+
+ // Decide whether down is positive
+ char posval[10];
+ int success = NCFUNC(get_att_text)(_fileId, (*vmit).second.varId, "positive", posval);
+ if (0 == success && !strcmp(posval, "down")) {
+ for (std::vector<double>::iterator dvit = levVals.begin(); dvit != levVals.end(); ++dvit)
+ (*dvit) *= -1.0;
+ }
+ }
+ else {
+ ERRORR(MB_FAILURE, "Couldn't find 'lev' variable.");
+ }
}
// Store time coordinate values in tVals
- if ((vmit = varInfo.find("time")) != varInfo.end() && (*vmit).second.varDims.size() == 1) {
- rval = read_coordinate("time", 0, nTimeSteps - 1, tVals);
- ERRORR(rval, "Trouble reading 'time' variable.");
- }
- else if ((vmit = varInfo.find("t")) != varInfo.end() && (*vmit).second.varDims.size() == 1) {
- rval = read_coordinate("t", 0, nTimeSteps - 1, tVals);
- ERRORR(rval, "Trouble reading 't' variable.");
+ if (isConnFile) {
+ // Connectivity file might not have time variable
}
else {
- // If expected time variable is not available, set dummy time coordinate values to tVals
- for (int t = 0; t < nTimeSteps; t++)
- tVals.push_back((double)t);
+ if ((vmit = varInfo.find("time")) != varInfo.end() && (*vmit).second.varDims.size() == 1) {
+ rval = read_coordinate("time", 0, nTimeSteps - 1, tVals);
+ ERRORR(rval, "Trouble reading 'time' variable.");
+ }
+ else if ((vmit = varInfo.find("t")) != varInfo.end() && (*vmit).second.varDims.size() == 1) {
+ rval = read_coordinate("t", 0, nTimeSteps - 1, tVals);
+ ERRORR(rval, "Trouble reading 't' variable.");
+ }
+ else {
+ // If expected time variable is not available, set dummy time coordinate values to tVals
+ for (int t = 0; t < nTimeSteps; t++)
+ tVals.push_back((double)t);
+ }
}
// For each variable, determine the entity location type and number of levels
@@ -225,23 +259,6 @@ ErrorCode NCHelperHOMME::create_mesh(Range& faces)
bool& spectralMesh = _readNC->spectralMesh;
int& gatherSetRank = _readNC->gatherSetRank;
- // Need to get/read connectivity data before creating elements
- std::string conn_fname;
-
- // Try to open the connectivity file through CONN option, if used
- ErrorCode rval = _opts.get_str_option("CONN", conn_fname);
- if (MB_SUCCESS != rval) {
- // Default convention for reading HOMME is a file HommeMapping.nc in same dir as data file
- conn_fname = std::string(fileName);
- size_t idx = conn_fname.find_last_of("/");
- if (idx != std::string::npos)
- conn_fname = conn_fname.substr(0, idx).append("/HommeMapping.nc");
- else
- conn_fname = "HommeMapping.nc";
- }
-
- int success = 0;
-
int rank = 0;
int procs = 1;
#ifdef USE_MPI
@@ -253,19 +270,42 @@ ErrorCode NCHelperHOMME::create_mesh(Range& faces)
}
#endif
+ ErrorCode rval;
+ int success = 0;
+
+ // Need to get/read connectivity data before creating elements
+ std::string conn_fname;
+
+ if (isConnFile) {
+ // Connectivity file has already been read
+ connectId = _readNC->fileId;
+ }
+ else {
+ // Try to open the connectivity file through CONN option, if used
+ rval = _opts.get_str_option("CONN", conn_fname);
+ if (MB_SUCCESS != rval) {
+ // Default convention for reading HOMME is a file HommeMapping.nc in same dir as data file
+ conn_fname = std::string(fileName);
+ size_t idx = conn_fname.find_last_of("/");
+ if (idx != std::string::npos)
+ conn_fname = conn_fname.substr(0, idx).append("/HommeMapping.nc");
+ else
+ conn_fname = "HommeMapping.nc";
+ }
#ifdef PNETCDF_FILE
#ifdef USE_MPI
- if (isParallel) {
- ParallelComm*& myPcomm = _readNC->myPcomm;
- success = NCFUNC(open)(myPcomm->proc_config().proc_comm(), conn_fname.c_str(), 0, MPI_INFO_NULL, &connectId);
- }
- else
- success = NCFUNC(open)(MPI_COMM_SELF, conn_fname.c_str(), 0, MPI_INFO_NULL, &connectId);
+ if (isParallel) {
+ ParallelComm*& myPcomm = _readNC->myPcomm;
+ success = NCFUNC(open)(myPcomm->proc_config().proc_comm(), conn_fname.c_str(), 0, MPI_INFO_NULL, &connectId);
+ }
+ else
+ success = NCFUNC(open)(MPI_COMM_SELF, conn_fname.c_str(), 0, MPI_INFO_NULL, &connectId);
#endif
#else
- success = NCFUNC(open)(conn_fname.c_str(), 0, &connectId);
+ success = NCFUNC(open)(conn_fname.c_str(), 0, &connectId);
#endif
- ERRORS(success, "Failed on open.");
+ ERRORS(success, "Failed on open.");
+ }
std::vector<std::string> conn_names;
std::vector<int> conn_vals;
@@ -286,7 +326,7 @@ ErrorCode NCHelperHOMME::create_mesh(Range& faces)
ERRORR(MB_FAILURE, "Failed to get number of quads.");
}
int num_quads = conn_vals[idx];
- if (num_quads != nCells) {
+ if (!isConnFile && num_quads != nCells) {
dbgOut.tprintf(1, "Warning: number of quads from %s and cells from %s are inconsistent; num_quads = %d, nCells = %d.\n",
conn_fname.c_str(), fileName.c_str(), num_quads, nCells);
}
@@ -300,8 +340,13 @@ ErrorCode NCHelperHOMME::create_mesh(Range& faces)
std::vector<int> tmp_conn(4 * num_quads), tmp_conn2(4 * num_quads);
success = NCFUNCAG(_vara_int)(connectId, cornerVarId, tmp_starts, tmp_counts, &tmp_conn2[0]);
ERRORS(success, "Failed to get temporary connectivity.");
- success = NCFUNC(close)(connectId);
- ERRORS(success, "Failed on close.");
+ if (isConnFile) {
+ // This data/connectivity file will be closed later in ReadNC::load_file()
+ }
+ else {
+ success = NCFUNC(close)(connectId);
+ ERRORS(success, "Failed on close.");
+ }
// Permute the connectivity
for (int i = 0; i < num_quads; i++) {
tmp_conn[4 * i] = tmp_conn2[i];
@@ -383,12 +428,12 @@ ErrorCode NCHelperHOMME::create_mesh(Range& faces)
// Convert lon/lat/rad to x/y/z
const double pideg = acos(-1.0) / 180.0;
+ double rad = (isConnFile) ? 8000.0 : 8000.0 + levVals[0];
for (i = 0; i < nLocalVertices; i++) {
double cosphi = cos(pideg * yptr[i]);
double zmult = sin(pideg * yptr[i]);
double xmult = cosphi * cos(xptr[i] * pideg);
double ymult = cosphi * sin(xptr[i] * pideg);
- double rad = 8000.0 + levVals[0];
xptr[i] = rad * xmult;
yptr[i] = rad * ymult;
zptr[i] = rad * zmult;
@@ -463,7 +508,6 @@ ErrorCode NCHelperHOMME::create_mesh(Range& faces)
double zmult = sin(pideg * yVertVals[i]);
double xmult = cosphi * cos(xVertVals[i] * pideg);
double ymult = cosphi * sin(xVertVals[i] * pideg);
- double rad = 8000.0 + levVals[0];
xptr[i] = rad * xmult;
yptr[i] = rad * ymult;
zptr[i] = rad * zmult;
diff --git a/src/io/NCHelperHOMME.hpp b/src/io/NCHelperHOMME.hpp
index 82ba90e..a6d50f0 100644
--- a/src/io/NCHelperHOMME.hpp
+++ b/src/io/NCHelperHOMME.hpp
@@ -46,6 +46,7 @@ private:
private:
int _spectralOrder; // Read from variable 'np'
int connectId; // For connectivity file
+ bool isConnFile; // Is the data file being read actually a connectivity file in disguise?
};
} // namespace moab
diff --git a/src/moab/SpatialLocator.hpp b/src/moab/SpatialLocator.hpp
index 0307b1c..f07cff7 100644
--- a/src/moab/SpatialLocator.hpp
+++ b/src/moab/SpatialLocator.hpp
@@ -133,6 +133,12 @@ namespace moab {
*/
Interface* moab() { return mbImpl; }
+ /* locate a point */
+ ErrorCode locate_point(const double *pos,
+ EntityHandle &ent, double *params, bool *is_inside = NULL,
+ const double rel_iter_tol = 1.0e-10, const double abs_iter_tol = 1.0e-10,
+ const double inside_tol = 1.0e-6);
+
/* return the tree */
Tree *get_tree() {return myTree;}
This diff is so big that we needed to truncate the remainder.
Repository URL: https://bitbucket.org/fathomteam/moab/
--
This is a commit notification from bitbucket.org. You are receiving
this because you have the service enabled, addressing the recipient of
this email.
More information about the moab-dev
mailing list