00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017
00018
00019
00020
00021
00022
00023
00024 static char const Id[] = "$Id: hstep_file.c 1000 2008-07-31 16:37:20Z tibus $\n\n";
00025 static char const Td[] = "$Today: " __FILE__ " " __DATE__ " " __TIME__ " $\n\n";
00026
00027 #include "field.h"
00028
00029 #define HSTEPCOL 2
00030
00031
00032
00033
00034
00035
00036
00037
00038
00039
00040
00041
00042
00043
00044
00045 int read_hstep_file(FILE *fd,PetscReal **indata,const int col)
00046 {
00047 int nrows;
00048 PetscReal *xydata;
00049
00050 MagparFunctionLogBegin;
00051
00052 int rank,size;
00053 ierr = MPI_Comm_size(PETSC_COMM_WORLD,&size);CHKERRQ(ierr);
00054 ierr = MPI_Comm_rank(PETSC_COMM_WORLD,&rank);CHKERRQ(ierr);
00055
00056 if (col<2) {
00057 SETERRQ1(PETSC_ERR_ARG_CORRUPT,"col=%i < 2!\n",col);
00058 }
00059
00060
00061 if (!rank) {
00062 char msg[256];
00063
00064
00065 fscanf(fd,"%i\n",&nrows);
00066 PetscPrintf(PETSC_COMM_WORLD,"Reading column %i, %i lines\n",col,nrows);
00067
00068
00069 nrows += 2;
00070
00071 ierr = PetscMalloc(HSTEPCOL*nrows*sizeof(PetscReal),&xydata);CHKERRQ(ierr);
00072
00073
00074 for (int i=1; i<nrows-1; i++) {
00075 int k;
00076
00077 #ifdef PETSC_USE_SINGLE
00078 k=fscanf(fd,"%f",xydata+i*HSTEPCOL+0);
00079 #else
00080 k=fscanf(fd,"%lf",xydata+i*HSTEPCOL+0);
00081 #endif
00082 if (k!=1) {
00083 SETERRQ2(PETSC_ERR_FILE_UNEXPECTED,
00084 "xydata in corrupt in row %i col %i\n",
00085 i,1
00086 );
00087 }
00088
00089
00090 for (int j=1; j<col; j++) {
00091 #ifdef PETSC_USE_SINGLE
00092 k=fscanf(fd,"%f",xydata+i*HSTEPCOL+1);
00093 #else
00094 k=fscanf(fd,"%lf",xydata+i*HSTEPCOL+1);
00095 #endif
00096
00097 if (k!=1) {
00098 SETERRQ2(PETSC_ERR_FILE_UNEXPECTED,
00099 "xydata in corrupt in row %i col %i\n",
00100 i,j
00101 );
00102 }
00103 }
00104
00105 fgets(msg,256,fd);
00106 }
00107
00108
00109 xydata[0]=PETSC_MIN;
00110 xydata[1]=0.0;
00111
00112
00113 xydata[HSTEPCOL*(nrows-1)+0]=PETSC_MAX;
00114 xydata[HSTEPCOL*(nrows-1)+1]=0.0;
00115
00116 PetscPrintf(PETSC_COMM_WORLD,
00117 "%14s %14s\n",
00118 "step","scale"
00119 );
00120
00121 for (int i=0; i<nrows; i++) {
00122 PetscPrintf(PETSC_COMM_WORLD,
00123 "%14e %14e\n",
00124 xydata[HSTEPCOL*i+0],
00125 xydata[HSTEPCOL*i+1]
00126 );
00127 }
00128 }
00129
00130
00131 ierr = MPI_Bcast(&nrows,1,MPI_INT,0,PETSC_COMM_WORLD);CHKERRQ(ierr);
00132 if (rank) {
00133 ierr = PetscMalloc(HSTEPCOL*nrows*sizeof(PetscReal),&xydata);CHKERRQ(ierr);
00134 }
00135 ierr = MPI_Bcast(xydata,nrows*HSTEPCOL,MPIU_SCALAR,0,PETSC_COMM_WORLD);CHKERRQ(ierr);
00136 *indata=xydata;
00137
00138 MagparFunctionLogReturn(0);
00139 }
00140
00141 PetscReal get_scalf_hstep(PetscReal xydata[],PetscReal xval)
00142 {
00143 int tint;
00144 PetscReal result;
00145
00146 tint=0;
00147 while (xval > xydata[HSTEPCOL*(tint+1)]) tint++;
00148
00149
00150 result=xydata[HSTEPCOL*tint+1];
00151
00152 return(result);
00153 }
00154