skip to primary navigationskip to content

Whittle Laboratory Intranet

Dept. of Engineering

Studying at Cambridge

C source code icon icem2padram3d_v1.c — C source code, 10 KB (10762 bytes)

File contents

#include <stdio.h>
#include <stdlib.h>
#include <string.h>

#define NBLK 500
#define NPCH 100
#define MAXLINE 1000

#define min(A,B) ( (A) < (B) ? (A):(B) )
#define max(A,B) ( (A) > (B) ? (A):(B) )

main()
{
    FILE *fid, *fptr;

    double *x, *y, *z;

    int count1, count2;
    int i, j, im, jm, km;
    int imax, jmax, kmax;
    int nblk, nblocks; 
    int ni[NBLK], nj[NBLK], nk[NBLK];
    int nb_nw, nb_nx;
    int ns_nw[3], ne_nw[3];
    int ns_nx[3], ne_nx[3];
    int is_nw, ie_nw, js_nw, je_nw, ks_nw, ke_nw;
    int is_nx, ie_nx, js_nx, je_nx, ks_nx, ke_nx;
    int bcflag, lblock1, lblock2;
    int l1lim, l2lim, m1lim1[2], m1lim2[2], n1lim1[2], n1lim2[2];

    char line[MAXLINE];
    char label[] = "\0";
    char ijk_nw[6], ijk_nx[6];
    char dir_nw[3], dir_nx[3];
    char sgn_nw[3], sgn_nx[3];
    char i_sgn_nw, j_sgn_nw, k_sgn_nw;
    char i_sgn_nx, j_sgn_nx, k_sgn_nx;
    char lface1, lface2, ldir1, ldir2, lspec1[2], lspec2[2]; 
    char bctype[] = "\0";
    char c_nw, c_nx;
    char grid_in[]= "\0";
    char grid_out[]="\0";
//
//  open the boundata file for writing
//  
    fptr = fopen("p3d.boundata","w");
    fclose(fptr);
//
//  open the topo file for reading    
//   
    printf("\n Opening info.topo ...\n");

    if ((fid = fopen("info.topo","r"))==NULL)
	    printf("File info.topo cannot be opened for reading\n\n");

    for (i = 0; i<2; i++)
	fgets(line,MAXLINE,fid);

    printf("\n Reading block size data ...\n\n");
			    
//
//  read in ni,nj,nk for each block until all blocks read
//  (i.e. until nblck=0)
//

    imax = 0;
    jmax = 0;
    kmax = 0;

    nblocks = 0;
    
//  if the string length is equal to 1, there is no more block data  
    while (strlen(fgets(line,MAXLINE,fid)) > 1){
    	sscanf(line,"%7s%d %*d %*d %*d %d %d %d",&label,&nblk,&im,&jm,&km);
	printf("Block number %5d ni,nj,nk = %5d%5d%5d\n",nblk,im,jm,km);

	nblocks++;
	ni[nblk-1] = im;
	nj[nblk-1] = jm;
	nk[nblk-1] = km;

	imax = max(im,imax);
	jmax = max(jm,jmax);
	kmax = max(km,kmax);
    }

    printf("\n Total number of blocks = %d\n",nblocks);

    if (nblocks > NBLK){
	    printf("\n Stopping because of too many blocks\n");
	    exit(-1);
    }

    printf("\n Reading connectivity / periodic data ...\n\n");
//
//  start of connectivity / periodic data reading for this block
//
    while (fgets(line,MAXLINE,fid) != NULL){
//
//     if an empty line is read exit while loop
//	    
       if (strlen(line) <= 1)
	   break;
//
//     if this is not connectivity / periodic data, all such data
//     has been read, therefore move on to processing
//
       if (strncmp(line,"# Connectivity",14) == 0)
	  bcflag = 0;   
       else if (strncmp(line,"# Periodic",10) == 0)
	  bcflag = 1;     
       else if (strncmp(line,"# Boundary",10) == 0){
	  bcflag = 2;
          sscanf(line,"%*51c%d",&nb_nw);
	  nb_nx = 0; //given blank value
       }	  
       else
	  break;

       printf("%s\n",line);
//
//     read in data for this patch
//     if the string length is equal to 1, there is no more 
//     connectivity information for this patch    
       while (fgets(line,MAXLINE,fid) != NULL){
	  if (strlen(line) <=1)
	     break;	  
	  if (bcflag < 2){     
             sscanf(line,"%*d %*7s%d %6c %c %d %d %d %d %d %d",&nb_nw,
			  &ijk_nw,&c_nw,&ns_nw[0],&ns_nw[1],&ns_nw[2],
			  &ne_nw[0],&ne_nw[1],&ne_nw[2]);
	     
//           before reading adjacent patch check for EOF
	     if (fgets(line,MAXLINE,fid) == NULL)
		 break;    
             sscanf(line,"%*d %*7s%d %6c %c %d %d %d %d %d %d",&nb_nx,
                          &ijk_nx,&c_nx,&ns_nx[0],&ns_nx[1],&ns_nx[2],
                          &ne_nx[0],&ne_nx[1],&ne_nx[2]);
	  }
	  if (bcflag == 2){
             sscanf(line,"%s %c %d %d %d %d %d %d",&bctype,&c_nw,&ns_nw[0]
			 ,&ns_nw[1],&ns_nw[2],&ne_nw[0],&ne_nw[1],&ne_nw[2]);
//
//           direction and sign are fixed
//	     
	     strcpy(ijk_nw," i j k");
	  }
//
//        scan through array ijk_nw and ijk_nx for i j and k
//	  
	  count1 = 0;
	  count2 = 0;
	  for (i = 0; i<6; i++){
	     if (ijk_nw[i]=='i' || ijk_nw[i]=='j' || ijk_nw[i]=='k'){
		     dir_nw[count1] = ijk_nw[i];
	             if (i == 0) sgn_nw[count1] = ' ';
		     else sgn_nw[count1] = ijk_nw[i-1];
	             count1++;
	     }
	     if (ijk_nx[i]=='i' || ijk_nx[i]=='j' || ijk_nx[i]=='k'){
		     dir_nx[count2] = ijk_nx[i];
		     if (i == 0) sgn_nx[count2] = ' ';
		     else sgn_nx[count2] = ijk_nx[i-1];
		     count2++;
	     }
	  }
//
//        if the current connectivity data is for a face (type f) then
//        padram needs this data for a patch - otherwise, skip to next
//        connectivity data
//
          if (c_nw == 'f'){
             for (i=0; i<3; i++){
		     if (dir_nw[i] == 'i'){
			 is_nw=min(ns_nw[i],ne_nw[i]);
		         ie_nw=max(ns_nw[i],ne_nw[i]);
			 i_sgn_nw = sgn_nw[i];}
		     else if (dir_nw[i] == 'j'){
			 js_nw=min(ns_nw[i],ne_nw[i]);    
		         je_nw=max(ns_nw[i],ne_nw[i]);
		         j_sgn_nw = sgn_nw[i];}
		     else if (dir_nw[i] == 'k'){
			 ks_nw=min(ns_nw[i],ne_nw[i]);
	                 ke_nw=max(ns_nw[i],ne_nw[i]);
		         k_sgn_nw = sgn_nw[i];}
//
                     if (dir_nx[i] == 'i'){
                         is_nx=min(ns_nx[i],ne_nx[i]);
                         ie_nx=max(ns_nx[i],ne_nx[i]);
		         i_sgn_nx = sgn_nx[i];}
                     else if (dir_nx[i] == 'j'){
                         js_nx=min(ns_nx[i],ne_nx[i]);
                         je_nx=max(ns_nx[i],ne_nx[i]);
		         j_sgn_nx = sgn_nx[i];}
                     else if (dir_nx[i] == 'k'){
                         ks_nx=min(ns_nx[i],ne_nx[i]);
                         ke_nx=max(ns_nx[i],ne_nx[i]);
			 k_sgn_nx = sgn_nx[i];}
	     }
//
//           connectivity data required by padram
//	     
	     lblock1 = nb_nw; lblock2 = nb_nx;
	     if (is_nw == ie_nw){
		 lface1 = 'I'; lspec1[0] = 'J'; lspec2[0] = 'K';
		 l1lim = is_nw; 
		 m1lim1[0] = js_nw; m1lim2[0] = je_nw;
		 n1lim1[0] = ks_nw; n1lim2[0] = ke_nw;
		 if (i_sgn_nw == ' ') ldir1 = 'P';
	         else ldir1 = 'M';
	     } 
	     else if (js_nw == je_nw){
		 lface1 = 'J'; lspec1[0] = 'I'; lspec2[0] = 'K';
		 l1lim = js_nw;
		 m1lim1[0] = is_nw; m1lim2[0] = ie_nw;
		 n1lim1[0] = ks_nw; n1lim2[0] = ke_nw;
		 if (j_sgn_nw == ' ') ldir1 = 'P';
	         else ldir1 = 'M';
	     }
	     else if (ks_nw == ke_nw){
		 lface1 = 'K'; lspec1[0] = 'I'; lspec2[0] = 'J';
		 l1lim = ks_nw;
		 m1lim1[0] = is_nw; m1lim2[0] = ie_nw;
		 n1lim1[0] = js_nw; n1lim2[0] = je_nw;
		 if (k_sgn_nw == ' ') ldir1 = 'P';
	         else ldir1 = 'M';
	     }
	     if (is_nx == ie_nx){
		 lface2 = 'I'; lspec1[1] = 'J'; lspec2[1] = 'K';
		 l2lim = is_nx;
		 m1lim1[1] = js_nx; m1lim2[1] = je_nx;
		 n1lim1[1] = ks_nx; n1lim2[1] = ke_nx;
		 if (i_sgn_nx == ' ') ldir2 = 'P';
	         else ldir2 = 'M';
	     }
	     else if (js_nx == je_nx){
		 lface2 = 'J'; lspec1[1] = 'I'; lspec2[1] = 'K';
		 l2lim = js_nx;
		 m1lim1[1] = is_nx; m1lim2[1] = ie_nx;
		 n1lim1[1] = ks_nx; n1lim2[1] = ke_nx;
		 if (j_sgn_nx == ' ') ldir2 = 'P';
		 else ldir2 = 'M';
	     }
	     else if (ks_nx == ke_nx){
		 lface2 = 'K'; lspec1[1] = 'I'; lspec2[1] = 'J'; 
		 l2lim = ks_nx;
		 m1lim1[1] = is_nx; m1lim2[1] = ie_nx;
		 n1lim1[1] = js_nx; n1lim2[1] = je_nx;
		 if (k_sgn_nx == ' ') ldir2 = 'P';
		 else ldir2 = 'M';
	     }
	    if (bcflag == 0)
		strcpy(bctype,"PATCH");
	    else if (bcflag == 1)
		strcpy(bctype,"PERIODIC"); 
//
//          write connectivity to p3d.boundata
//	  
	    fptr = fopen("p3d.boundata","a");
	    if (bcflag != 2)
		if (nb_nw != nb_nx)    
                    fprintf(fptr,"%-10s%4d%4d%3c%3c%3c%3c%3c%3c%4d%4d%4d%4d"
			 "%4d%4d%4d%4d%4d%4d\n",bctype,lblock1,lblock2
			 ,lface1,lface2,ldir1,ldir2,lspec1[0],lspec2[0],
                         l1lim,l2lim,m1lim1[0],m1lim2[0],n1lim1[0],n1lim2[0],
			 m1lim1[0],m1lim2[0],n1lim1[0],n1lim2[0]);  
	        else
		{
		    fprintf(fptr,"%-10s%4d%4d%3c%3c%3c%3c%3c%3c%4d%4d%4d%4d"
			 "%4d%4d%4d%4d%4d%4d\n",bctype,lblock1,lblock2
		         ,lface1,lface1,ldir1,ldir1,lspec1[0],lspec2[0],
		         l1lim,l1lim,m1lim1[0],m1lim2[0],n1lim1[0],n1lim2[0],
		         m1lim1[0],m1lim2[0],n1lim1[0],n1lim2[0]);
                    fprintf(fptr,"%-10s%4d%4d%3c%3c%3c%3c%3c%3c%4d%4d%4d%4d"
		         "%4d%4d%4d%4d%4d%4d\n",bctype,lblock1,lblock2
		         ,lface2,lface2,ldir2,ldir2,lspec1[1],lspec2[1],
		         l2lim,l2lim,m1lim1[1],m1lim2[1],n1lim1[1],n1lim2[1],
		         m1lim1[1],m1lim2[1],n1lim1[1],n1lim2[1]);	 
		}
	    else if (bcflag == 2)
	        fprintf(fptr,"%-10s%4d%4d%3c%3c%3c%3c%3c%3c%4d%4d%4d%4d"
			 "%4d%4d%4d%4d%4d%4d\n",bctype,lblock1,lblock1
			 ,lface1,lface1,ldir1,ldir1,lspec1[0],lspec2[0],
			 l1lim,l1lim,m1lim1[0],m1lim2[0],n1lim1[0],n1lim2[0],
			 m1lim1[0],m1lim2[0],n1lim1[0],n1lim2[0]);
//	    
	    fclose(fptr);
	  }		  
       }
    }

    printf("\n Closing info.topo ...\n");
    fclose(fid);
    
//
//  all connectivity / periodicity data has now been read in
//    

    printf("\n Creating mesh.plot3d from info.dom files ...\n");
//    
    sprintf(grid_out,"%s","mesh.plot3d");
    fid = fopen(grid_out,"w");
//
//  write header line for mesh.plot3d file
//  writing number of blocks
//
    fprintf(fid,"%4d\n",nblocks);

//  write block dimensions
    for (i=0; i<nblocks; i++)
	fprintf(fid,"%5d %5d %5d",ni[i],nj[i],nk[i]);  
    fprintf(fid,"\n");  
	    
//
//  read info.dom files
//    
    for (i=0; i<nblocks; i++){ 
    	sprintf(grid_in,"%s%d","info.dom",i+1);
	printf("Opening file %s\n",grid_in);

	if ((fptr = fopen(grid_in,"r")) == NULL)
		return -1;
//      read headerline
	fgets(line,MAXLINE,fptr);
//
//      allocate memory for node locations
//	
	x = (double *) malloc(sizeof(double) * (ni[i]*nj[i]*nk[i]));
	y = (double *) malloc(sizeof(double) * (ni[i]*nj[i]*nk[i]));
	z = (double *) malloc(sizeof(double) * (ni[i]*nj[i]*nk[i]));
	
	for (j=0; j<ni[i]*nj[i]*nk[i]; j++)
	    fscanf(fptr,"%lf %lf %lf\n",&x[j],&y[j],&z[j]);
	fclose(fptr);
//	    
//      write out all the x followed by all the y and z coords
//	    
        for (j=0; j<ni[i]*nj[i]*nk[i]; j++)
	    fprintf(fid,"%16.8f",x[j]);
 	for (j=0; j<ni[i]*nj[i]*nk[i]; j++)
	    fprintf(fid,"%16.8f",y[j]);
	for (j=0; j<ni[i]*nj[i]*nk[i]; j++)
	    fprintf(fid,"%16.8f",z[j]);
	fprintf(fid,"\n");
//
//      free memory
//	
	free(x);
	free(y);
	free(z);
    }
    
    fclose(fid);
}
« July 2017 »
July
MoTuWeThFrSaSu
12
3456789
10111213141516
17181920212223
24252627282930
31