Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Correction of reading .npy files #13

Merged
merged 2 commits into from
Feb 16, 2024
Merged

Correction of reading .npy files #13

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
2 commits
Select commit Hold shift + click to select a range
08015b4
Added output messages in npy reading
Jan 31, 2024
e5f6edf
Proper input of initial data from .npy file : number of entities in e…
Feb 1, 2024
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
2 changes: 1 addition & 1 deletion src/hashhexa.c
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -249,7 +249,7 @@ int H2T_hashQuad(MMG5_pMesh mesh,MMG5_Hash *hash) {
*/
int H2T_hashGetRef(MMG5_pMesh mesh,MMG5_Hash *hash) {
MMG5_pTetra pt;
int *adja;
MMG5_int *adja;
int ie,i,nt,ia,ib,ic,k,ref;

/* Create tetrahedra adjacency */
Expand Down
86 changes: 46 additions & 40 deletions src/inout.c
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -16,17 +16,17 @@
#include <math.h>

static inline
int H2T_npy_point_index(int i,int j,int k, int *n) {
return n[1]*n[2]*i + n[2]*j + k + 1;
MMG5_int H2T_npy_point_index(int i,int j,int k, int *n) {
return (n[1]+1)*(n[2]+1)*i + (n[2]+1)*j + k + 1;
}

int H2T_loadNpy(MMG5_pMesh mmgMesh, int** tabhex, char* filename) {

FILE* inm;
unsigned char buffer = 0x00;
char* str = NULL;
int pos1, pos2, dim = 0, t[3], nhex;
MMG5_int np, ne, i, j, k, ref, pos;
int pos1, pos2, dim = 0, t[3];
MMG5_int nhex, np, ne, i, j, k, ref, pos;

/* Input data and creation of the hexa array */
if( !(inm = fopen(mmgMesh->namein,"rb")) ) {
Expand Down Expand Up @@ -71,9 +71,9 @@ int H2T_loadNpy(MMG5_pMesh mmgMesh, int** tabhex, char* filename) {
fread(&buffer,sizeof buffer,1,inm);
}

np = t[0]*t[1]*t[2];
nhex = (t[0]-1)*(t[1]-1)*(t[2]-1);
ne = 4 * ((t[0]-1)+(t[1]-1)+(t[2]-1));
np = (t[0]+1)*(t[1]+1)*(t[2]+1);
nhex = t[0]*t[1]*t[2];
ne = 4 * (t[0]+t[1]+t[2]);

if ( H2T_Set_meshSize(mmgMesh,np,nhex,0,ne) != 1 ) {
return -1;
Expand All @@ -84,9 +84,10 @@ int H2T_loadNpy(MMG5_pMesh mmgMesh, int** tabhex, char* filename) {
/* Point coordinates from grid indices */
ref = 0;
pos = 0;
for (i=0;i<t[0];i++) {
for (j=0;j<t[1];j++) {
for (k=0;k<t[2];k++) {
fprintf(stdout," READING %" MMG5_PRId " VERTICES\n",np);
for (i=0;i<=t[0];i++) {
for (j=0;j<=t[1];j++) {
for (k=0;k<=t[2];k++) {
if ( H2T_Set_vertex(mmgMesh,(double)i ,(double)j ,(double)k,ref,++pos) != 1 ) {
return -1;
}
Expand All @@ -96,10 +97,11 @@ int H2T_loadNpy(MMG5_pMesh mmgMesh, int** tabhex, char* filename) {

/* Hexahedra */
pos = 1;
for ( i=0; i<t[0]-1; ++i ) {
for ( j=0; j<t[1]-1; ++j ) {
for ( k=0; k<t[2]-1; ++k ) {
int iadr = 9*pos;
fprintf(stdout," READING %" MMG5_PRId " HEXA\n",nhex);
for ( i=0; i<t[0]; ++i ) {
for ( j=0; j<t[1]; ++j ) {
for ( k=0; k<t[2]; ++k ) {
MMG5_int iadr = 9*pos;

/* Hexa vertices */
(*tabhex)[iadr+0] = H2T_npy_point_index(i ,j ,k ,t);
Expand All @@ -112,7 +114,7 @@ int H2T_loadNpy(MMG5_pMesh mmgMesh, int** tabhex, char* filename) {
(*tabhex)[iadr+7] = H2T_npy_point_index(i ,j+1,k+1,t);

/* Hexa references */
fread(&(*tabhex)[iadr+8],sizeof(int16_t),1,inm);
fread(&(*tabhex)[iadr+8],sizeof(uint32_t),1,inm);
++pos;
}
}
Expand All @@ -121,64 +123,68 @@ int H2T_loadNpy(MMG5_pMesh mmgMesh, int** tabhex, char* filename) {
/* Edges */
ref = 0;
pos = 0;
for (i=0;i<t[0]-1; ++i) {
int np0, np1;
fprintf(stdout," READING %" MMG5_PRId " EDGES\n",ne);
for (i=0;i<t[0]; ++i) {
MMG5_int np0, np1;
np0 = H2T_npy_point_index(i ,0 ,0 ,t);
np1 = H2T_npy_point_index(i+1,0 ,0 ,t);
H2T_Set_edge(mmgMesh,np0,np1,ref,++pos);
if (H2T_Set_edge(mmgMesh,np0,np1,ref,++pos) != 1) {
return -1;
}

np0 = H2T_npy_point_index(i ,t[1]-1,0 ,t);
np1 = H2T_npy_point_index(i+1,t[1]-1,0 ,t);
np0 = H2T_npy_point_index(i ,t[1],0 ,t);
np1 = H2T_npy_point_index(i+1,t[1],0 ,t);
H2T_Set_edge(mmgMesh,np0,np1,ref,++pos);

np0 = H2T_npy_point_index(i ,0 ,t[2]-1,t);
np1 = H2T_npy_point_index(i+1,0 ,t[2]-1,t);
np0 = H2T_npy_point_index(i ,0 ,t[2],t);
np1 = H2T_npy_point_index(i+1,0 ,t[2],t);
H2T_Set_edge(mmgMesh,np0,np1,ref,++pos);

np0 = H2T_npy_point_index(i ,t[1]-1,t[2]-1,t);
np1 = H2T_npy_point_index(i+1,t[1]-1,t[2]-1,t);
np0 = H2T_npy_point_index(i ,t[1],t[2],t);
np1 = H2T_npy_point_index(i+1,t[1],t[2],t);
H2T_Set_edge(mmgMesh,np0,np1,ref,++pos);
}

for (j=0;j<t[1]-1; ++j) {
int np0, np1;
for (j=0;j<t[1]; ++j) {
MMG5_int np0, np1;
np0 = H2T_npy_point_index(0 ,j ,0 ,t);
np1 = H2T_npy_point_index(0 ,j+1,0 ,t);
H2T_Set_edge(mmgMesh,np0,np1,ref,++pos);

np0 = H2T_npy_point_index(t[0]-1,j ,0 ,t);
np1 = H2T_npy_point_index(t[0]-1,j+1,0 ,t);
np0 = H2T_npy_point_index(t[0],j ,0 ,t);
np1 = H2T_npy_point_index(t[0],j+1,0 ,t);
H2T_Set_edge(mmgMesh,np0,np1,ref,++pos);

np0 = H2T_npy_point_index(0 ,j ,t[2]-1,t);
np1 = H2T_npy_point_index(0 ,j+1,t[2]-1,t);
np0 = H2T_npy_point_index(0 ,j ,t[2],t);
np1 = H2T_npy_point_index(0 ,j+1,t[2],t);
H2T_Set_edge(mmgMesh,np0,np1,ref,++pos);

np0 = H2T_npy_point_index(t[0]-1,j ,t[2]-1,t);
np1 = H2T_npy_point_index(t[0]-1,j+1,t[2]-1,t);
np0 = H2T_npy_point_index(t[0],j ,t[2],t);
np1 = H2T_npy_point_index(t[0],j+1,t[2],t);
H2T_Set_edge(mmgMesh,np0,np1,ref,++pos);
}

for (k=0;k<t[2]-1; ++k) {
int np0, np1;
for (k=0;k<t[2]; ++k) {
MMG5_int np0, np1;
np0 = H2T_npy_point_index(0 ,0 ,k ,t);
np1 = H2T_npy_point_index(0 ,0 ,k+1,t);
H2T_Set_edge(mmgMesh,np0,np1,ref,++pos);

np0 = H2T_npy_point_index(t[0]-1,0 ,k ,t);
np1 = H2T_npy_point_index(t[0]-1,0 ,k+1,t);
np0 = H2T_npy_point_index(t[0],0 ,k ,t);
np1 = H2T_npy_point_index(t[0],0 ,k+1,t);
H2T_Set_edge(mmgMesh,np0,np1,ref,++pos);

np0 = H2T_npy_point_index(0 ,t[1]-1,k ,t);
np1 = H2T_npy_point_index(0 ,t[1]-1,k+1,t);
np0 = H2T_npy_point_index(0 ,t[1],k ,t);
np1 = H2T_npy_point_index(0 ,t[1],k+1,t);
H2T_Set_edge(mmgMesh,np0,np1,ref,++pos);

np0 = H2T_npy_point_index(t[0]-1,t[1]-1,k ,t);
np1 = H2T_npy_point_index(t[0]-1,t[1]-1,k+1,t);
np0 = H2T_npy_point_index(t[0],t[1],k ,t);
np1 = H2T_npy_point_index(t[0],t[1],k+1,t);
H2T_Set_edge(mmgMesh,np0,np1,ref,++pos);
}

/* Ridges */
fprintf(stdout," READING %" MMG5_PRId " RIDGES\n",ne);
for (i=1;i<=ne;i++) {
MMG3D_Set_ridge(mmgMesh,i);
}
Expand Down
Loading