HyPar  1.0
Finite-Difference Hyperbolic-Parabolic PDE Solver on Cartesian Grids
tridiagLUInit.c File Reference

Initialize the tridiagLU solver. More...

#include <stdio.h>
#include <string.h>
#include <tridiagLU.h>
#include <mpi.h>

Go to the source code of this file.

Functions

int tridiagLUInit (void *r, void *c)
 

Detailed Description

Initialize the tridiagLU solver.

Author
Debojyoti Ghosh

Definition in file tridiagLUInit.c.

Function Documentation

◆ tridiagLUInit()

int tridiagLUInit ( void *  r,
void *  c 
)

Initialize the tridiagLU solver by setting the various parameters in TridiagLU to their default values. If the file lusolver.inp is available, read it and set the parameters.

The file lusolver.inp must be in the following format:

    begin
        <keyword>   <value>
        <keyword>   <value>
        <keyword>   <value>
        ...
        <keyword>   <value>
    end

where the list of keywords are:

Keyword name Type Variable Default value ---—
evaluate_norm int TridiagLU::evaluate_norm 1
maxiter int TridiagLU::maxiter 10
atol double TridiagLU::atol 1e-12
rtol double TridiagLU::rtol 1e-10
verbose int TridiagLU::verbose 0
reducedsolvetype char[] TridiagLU::reducedsolvetype _TRIDIAG_JACOBI_
Parameters
rObject of type TridiagLU
cMPI communicator

Definition at line 39 of file tridiagLUInit.c.

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 }
char reducedsolvetype[50]
Definition: tridiagLU.h:91
int verbose
Definition: tridiagLU.h:99
int evaluate_norm
Definition: tridiagLU.h:93
double atol
Definition: tridiagLU.h:95
#define _TRIDIAG_JACOBI_
Definition: tridiagLU.h:74
double rtol
Definition: tridiagLU.h:95
int maxiter
Definition: tridiagLU.h:94