Timothy  0.9
Tissue Modelling Framework
 All Data Structures Files Functions Variables Typedefs Macros
main.c
Go to the documentation of this file.
1 /* **************************************************************************
2  * This file is part of Timothy
3  *
4  * Copyright (c) 2014/15 Maciej Cytowski
5  * Copyright (c) 2014/15 ICM, University of Warsaw, Poland
6  *
7  * This program is free software; you can redistribute it and/or modify
8  * it under the terms of the GNU General Public License as published by
9  * the Free Software Foundation; either version 2 of the License, or
10  * (at your option) any later version.
11  *
12  * This program is distributed in the hope that it will be useful,
13  * but WITHOUT ANY WARRANTY; without even the implied warranty of
14  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
15  * GNU General Public License for more details.
16  *
17  * You should have received a copy of the GNU General Public License
18  * along with this program; if not, write to the Free Software
19  * Foundation, Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
20  *
21  * *************************************************************************/
22 
23 #include<stdio.h>
24 #include<stdlib.h>
25 #include<string.h>
26 #include<mpi.h>
27 #include<omp.h>
28 
29 #include "global.h"
30 
40 int main(int argc, char **argv)
41 {
42 
43  MPI_Init(&argc, &argv);
44  MPI_Comm_size(MPI_COMM_WORLD, &MPIsize);
45  MPI_Comm_rank(MPI_COMM_WORLD, &MPIrank);
46 
47  OMPthreads = omp_get_max_threads();
48 
49  simulationInit(argc, argv);
50 
51  for (step = 0; step < nsteps; step++) {
52 
53  if (!(step % statOutStep))
54  printStepNum();
55 
57  treeBuild();
62 
64  fieldsSolve();
66 
67  if (!(step % statOutStep))
69 
70  if (simStart)
71  simTime += secondsPerStep / 3600.0; /* biological process time in hours */
72 
73  if (!(step % vtkOutStep)) {
74  if (vtkout)
76  if (povout)
78  if (vnfout)
80  }
81 
84  commCleanup();
85  treeFree();
86 
87  if (!(step % rstOutStep))
88  saveRstFile();
89  }
90 
91  MPI_Barrier(MPI_COMM_WORLD);
92 
95  cellsCleanup();
96 
97  if (MPIrank == 0)
98  printf("\nEnd of simulation run.\n");
99 
100  MPI_Finalize();
101 
102  return 0;
103 }
int statOutStep
Definition: global.h:214
void ioWriteStepPovRay(int step, int type)
Definition: io.c:1978
void ioWriteFields(int step)
Definition: io.c:1672
float secondsPerStep
Definition: global.h:247
void initCellsToGridExchange()
Definition: interp.c:598
int MPIrank
Definition: global.h:134
int povout
Definition: global.h:181
void treeBuild()
Definition: tree.c:181
int simStart
Definition: global.h:172
void treeFree()
Definition: tree.c:369
void decompositionFinalize()
Definition: domdec.c:277
int vtkOutStep
Definition: global.h:216
int MPIsize
Definition: global.h:135
void commCleanup()
Definition: comm.c:156
int nsteps
Definition: global.h:169
void waitCellsToGridExchange()
Definition: interp.c:613
void printStepNum()
Definition: io.c:1215
void saveRstFile()
Definition: io.c:1358
contains the most important global variables, arrays and defines
void randomStreamFree()
Definition: random.c:45
void decompositionExecute()
Definition: domdec.c:240
void updateCellPositions()
Definition: cells.c:504
void createExportList()
Definition: comm.c:64
void computePotential()
Definition: potential.c:287
int rstOutStep
Definition: global.h:215
int vtkout
Definition: global.h:180
int OMPthreads
Definition: global.h:137
int step
Definition: global.h:173
void cellsCleanup()
Definition: cells.c:434
void fieldsSolve()
Definition: fields.c:131
void statisticsPrint()
Definition: stats.c:153
float simTime
Definition: global.h:175
int vnfout
Definition: global.h:182
void updateCellStates()
Definition: cells.c:718
void ioWriteStepVTK(int step)
Definition: io.c:1001
void simulationInit(int argc, char **argv)
Definition: init.c:74
void interpolateFieldsToCells()
Definition: interp.c:627
int main(int argc, char **argv)
Definition: main.c:40
int computePotentialGradient()
Definition: gradient.c:183