HyPar  1.0
Finite-Difference Hyperbolic-Parabolic PDE Solver on Cartesian Grids
 All Data Structures Namespaces Files Functions Variables Typedefs Enumerations Enumerator Macros Pages
tridiagLUInit.c
Go to the documentation of this file.
1 
6 #include <stdio.h>
7 #include <string.h>
8 #include <tridiagLU.h>
9 #ifndef serial
10 #include <mpi.h>
11 #endif
12 
40  void *r,
41  void *c
42  )
43 {
44  TridiagLU *t = (TridiagLU*) r;
45  int rank,ierr;
46 #ifdef serial
47  rank = 0;
48 #else
49  MPI_Comm *comm = (MPI_Comm*) c;
50  if (!comm) rank = 0;
51  else MPI_Comm_rank(*comm,&rank);
52 #endif
53 
54  /* default values */
56  t->evaluate_norm = 1;
57  t->maxiter = 10;
58  t->atol = 1e-12;
59  t->rtol = 1e-10;
60  t->verbose = 0;
61 
62  /* read from file, if available */
63  if (!rank) {
64  FILE *in;
65  in = fopen("lusolver.inp","r");
66  if (!in) {
67  printf("tridiagLUInit: File \"lusolver.inp\" not found. Using default values.\n");
68  } else {
69  char word[100];
70  ierr = fscanf(in,"%s",word); if (ierr != 1) return(1);
71  if (!strcmp(word, "begin")) {
72  while (strcmp(word, "end")) {
73  ierr = fscanf(in,"%s",word); if (ierr != 1) return(1);
74  if (!strcmp(word, "evaluate_norm" )) ierr = fscanf(in,"%d" ,&t->evaluate_norm );
75  else if (!strcmp(word, "maxiter" )) ierr = fscanf(in,"%d" ,&t->maxiter );
76  else if (!strcmp(word, "atol" )) ierr = fscanf(in,"%lf",&t->atol );
77  else if (!strcmp(word, "rtol" )) ierr = fscanf(in,"%lf",&t->rtol );
78  else if (!strcmp(word, "verbose" )) ierr = fscanf(in,"%d" ,&t->verbose );
79  else if (!strcmp(word, "reducedsolvetype")) ierr = fscanf(in,"%s" ,t->reducedsolvetype);
80  else if (strcmp(word,"end")) {
81  char useless[100];
82  ierr = fscanf(in,"%s",useless);
83  printf("Warning: keyword %s in file \"lusolver.inp\" with value %s not recognized or extraneous. Ignoring.\n",word,useless);
84  }
85  if (ierr != 1) return(1);
86  }
87  } else {
88  fprintf(stderr,"Error: Illegal format in file \"solver.inp\".\n");
89  return(1);
90  }
91  fclose(in);
92  }
93  }
94 
95  /* broadcast to all processes */
96 #ifndef serial
97  if (comm) {
98  MPI_Bcast(t->reducedsolvetype,50,MPI_CHAR,0,*comm);
99  MPI_Bcast(&t->evaluate_norm,1,MPI_INT,0,*comm);
100  MPI_Bcast(&t->maxiter,1,MPI_INT,0,*comm);
101  MPI_Bcast(&t->verbose,1,MPI_INT,0,*comm);
102  MPI_Bcast(&t->atol,1,MPI_DOUBLE,0,*comm);
103  MPI_Bcast(&t->rtol,1,MPI_DOUBLE,0,*comm);
104  }
105 #endif
106 
107  return(0);
108 }
109 
Header file for TridiagLU.
int verbose
Definition: tridiagLU.h:99
double rtol
Definition: tridiagLU.h:95
char reducedsolvetype[50]
Definition: tridiagLU.h:91
double atol
Definition: tridiagLU.h:95
int tridiagLUInit(void *, void *)
Definition: tridiagLUInit.c:39
int maxiter
Definition: tridiagLU.h:94
int evaluate_norm
Definition: tridiagLU.h:93
#define _TRIDIAG_JACOBI_
Definition: tridiagLU.h:74