| 1 |
/* $Id:$ */ |
| 2 |
|
| 3 |
/******************************************************* |
| 4 |
* |
| 5 |
* Copyright 2007 by University of Queensland |
| 6 |
* |
| 7 |
* http://esscc.uq.edu.au |
| 8 |
* Primary Business: Queensland, Australia |
| 9 |
* Licensed under the Open Software License version 3.0 |
| 10 |
* http://www.opensource.org/licenses/osl-3.0.php |
| 11 |
* |
| 12 |
*******************************************************/ |
| 13 |
|
| 14 |
/**************************************************************/ |
| 15 |
|
| 16 |
/* Paso: FCTransportProblem */ |
| 17 |
|
| 18 |
/**************************************************************/ |
| 19 |
|
| 20 |
/* Author: l.gross@uq.edu.au */ |
| 21 |
|
| 22 |
/**************************************************************/ |
| 23 |
|
| 24 |
|
| 25 |
#include "Paso.h" |
| 26 |
#include "SolverFCT.h" |
| 27 |
#include "PasoUtil.h" |
| 28 |
|
| 29 |
|
| 30 |
/**************************************************************/ |
| 31 |
|
| 32 |
/* free all memory used by */ |
| 33 |
void Paso_FCTransportProblem_free(Paso_FCTransportProblem* in) { |
| 34 |
if (in!=NULL) { |
| 35 |
in->reference_counter--; |
| 36 |
if (in->reference_counter<=0) { |
| 37 |
|
| 38 |
Paso_SystemMatrix_free(in->transport_matrix); |
| 39 |
Paso_SystemMatrix_free(in->flux_matrix); |
| 40 |
Paso_MPIInfo_free(in->mpi_info); |
| 41 |
|
| 42 |
MEMFREE(in->u); |
| 43 |
MEMFREE(in->lumped_mass_matrix); |
| 44 |
MEMFREE(in->row_sum_flux_matrix); |
| 45 |
MEMFREE(in->colorOf); |
| 46 |
MEMFREE(in->main_iptr); |
| 47 |
MEMFREE(in); |
| 48 |
} |
| 49 |
} |
| 50 |
} |
| 51 |
|
| 52 |
Paso_FCTransportProblem* Paso_FCTransportProblem_getReference(Paso_FCTransportProblem* in) { |
| 53 |
if (in!=NULL) { |
| 54 |
++(in->reference_counter); |
| 55 |
} |
| 56 |
return in; |
| 57 |
} |
| 58 |
|
| 59 |
Paso_SystemMatrix* Paso_FCTransportProblem_borrowTransportMatrix(Paso_FCTransportProblem* in) { |
| 60 |
return in->transport_matrix; |
| 61 |
} |
| 62 |
|
| 63 |
Paso_SystemMatrix* Paso_FCTransportProblem_borrowFluxMatrix(Paso_FCTransportProblem* in) { |
| 64 |
return in->flux_matrix; |
| 65 |
} |
| 66 |
|
| 67 |
double* Paso_FCTransportProblem_borrowLumpedMassMatrix(Paso_FCTransportProblem* in) { |
| 68 |
return in->lumped_mass_matrix; |
| 69 |
} |
| 70 |
|
| 71 |
dim_t Paso_FCTransportProblem_getTotalNumRows(Paso_FCTransportProblem* in) { |
| 72 |
return Paso_SystemMatrix_getTotalNumRows(in->transport_matrix); |
| 73 |
} |
| 74 |
|
| 75 |
Paso_FCTransportProblem* Paso_FCTransportProblem_alloc(double theta, double dt_max, Paso_SystemMatrixPattern *pattern, int block_size |
| 76 |
|
| 77 |
|
| 78 |
) { |
| 79 |
Paso_SystemMatrixType matrix_type=MATRIX_FORMAT_DEFAULT+MATRIX_FORMAT_BLK1; /* at the moment only block size 1 is supported */ |
| 80 |
Paso_FCTransportProblem* out=NULL; |
| 81 |
dim_t n,i; |
| 82 |
index_t iptr,iptr_main,k; |
| 83 |
|
| 84 |
if ((theta<0.) || (theta >1.)) { |
| 85 |
Paso_setError(TYPE_ERROR,"Paso_FCTransportProblem_alloc: theta needs to be between 0. and. 1."); |
| 86 |
return NULL; |
| 87 |
} |
| 88 |
|
| 89 |
out=MEMALLOC(1,Paso_FCTransportProblem); |
| 90 |
if (Paso_checkPtr(out)) return NULL; |
| 91 |
|
| 92 |
out->theta=theta; |
| 93 |
out->valid_matrices=FALSE; |
| 94 |
out->transport_matrix=Paso_SystemMatrix_alloc(matrix_type,pattern,block_size,block_size); |
| 95 |
Paso_SystemMatrix_allocBuffer(out->transport_matrix); |
| 96 |
out->flux_matrix=Paso_SystemMatrix_alloc(matrix_type,pattern,block_size,block_size); |
| 97 |
out->mpi_info=Paso_MPIInfo_getReference(pattern->mpi_info); |
| 98 |
|
| 99 |
out->colorOf=NULL; |
| 100 |
out->main_iptr=NULL; |
| 101 |
out->lumped_mass_matrix=NULL; |
| 102 |
out->row_sum_flux_matrix=NULL; |
| 103 |
|
| 104 |
if (Paso_noError()) { |
| 105 |
n=Paso_SystemMatrix_getTotalNumRows(out->transport_matrix); |
| 106 |
|
| 107 |
out->colorOf=MEMALLOC(n,index_t); |
| 108 |
out->main_iptr=MEMALLOC(n,index_t); |
| 109 |
out->lumped_mass_matrix=MEMALLOC(n,double); |
| 110 |
out->row_sum_flux_matrix=MEMALLOC(n,double); |
| 111 |
out->u=MEMALLOC(n,double); |
| 112 |
|
| 113 |
if ( ! (Paso_checkPtr(out->colorOf) || Paso_checkPtr(out->main_iptr) || |
| 114 |
Paso_checkPtr(out->lumped_mass_matrix) || Paso_checkPtr(out->row_sum_flux_matrix) || Paso_checkPtr(out->u)) ) { |
| 115 |
|
| 116 |
printf("Paso_SolverFCT_getFCTransportProblem: Revise coloring!!\n"); |
| 117 |
Paso_Pattern_color(pattern->mainPattern,&(out->num_colors),out->colorOf); |
| 118 |
|
| 119 |
|
| 120 |
#pragma omp parallel for schedule(static) private(i) |
| 121 |
for (i = 0; i < n; ++i) { |
| 122 |
out->lumped_mass_matrix[i]=0.; |
| 123 |
out->row_sum_flux_matrix[i]=0.; |
| 124 |
out->u[i]=0.; |
| 125 |
} |
| 126 |
|
| 127 |
/* identify the main diagonals */ |
| 128 |
#pragma omp parallel for schedule(static) private(i,iptr,iptr_main,k) |
| 129 |
for (i = 0; i < n; ++i) { |
| 130 |
for (iptr=pattern->mainPattern->ptr[i];iptr<pattern->mainPattern->ptr[i+1]; ++iptr) { |
| 131 |
iptr_main=pattern->mainPattern->ptr[0]-1; |
| 132 |
for (iptr=pattern->mainPattern->ptr[i];iptr<pattern->mainPattern->ptr[i+1]; iptr++) { |
| 133 |
if (pattern->mainPattern->index[iptr]==i) { |
| 134 |
iptr_main=iptr; |
| 135 |
break; |
| 136 |
} |
| 137 |
} |
| 138 |
out->main_iptr[i]=iptr_main; |
| 139 |
if (iptr_main==pattern->mainPattern->ptr[0]-1) |
| 140 |
Paso_setError(VALUE_ERROR, "Paso_FCTransportProblem_alloc: no main diagonal"); |
| 141 |
} |
| 142 |
} |
| 143 |
|
| 144 |
} |
| 145 |
|
| 146 |
} |
| 147 |
if (Paso_noError()) { |
| 148 |
return out; |
| 149 |
} else { |
| 150 |
Paso_FCTransportProblem_free(out); |
| 151 |
return NULL; |
| 152 |
} |
| 153 |
} |
| 154 |
|
| 155 |
void Paso_FCTransportProblem_checkinSolution(Paso_FCTransportProblem* in, double* u) { |
| 156 |
dim_t i, n; |
| 157 |
|
| 158 |
n=Paso_FCTransportProblem_getTotalNumRows(in); |
| 159 |
#pragma omp parallel for schedule(static) private(i) |
| 160 |
for (i = 0; i < n; ++i) { |
| 161 |
in->u[i]=u[i]; |
| 162 |
} |
| 163 |
} |
Parent Directory
| 
Revision Log