/[escript]/branches/doubleplusgood/paso/src/PCG.cpp
ViewVC logotype

Diff of /branches/doubleplusgood/paso/src/PCG.cpp

Parent Directory Parent Directory | Revision Log Revision Log | View Patch Patch

revision 1579 by ksteube, Mon Jun 2 08:48:36 2008 UTC revision 1811 by ksteube, Thu Sep 25 23:11:13 2008 UTC
# Line 1  Line 1 
1    
 /* $Id$ */  
   
2  /*******************************************************  /*******************************************************
3   *  *
4   *           Copyright 2003-2007 by ACceSS MNRF  * Copyright (c) 2003-2008 by University of Queensland
5   *       Copyright 2007 by University of Queensland  * Earth Systems Science Computational Center (ESSCC)
6   *  * http://www.uq.edu.au/esscc
7   *                http://esscc.uq.edu.au  *
8   *        Primary Business: Queensland, Australia  * Primary Business: Queensland, Australia
9   *  Licensed under the Open Software License version 3.0  * Licensed under the Open Software License version 3.0
10   *     http://www.opensource.org/licenses/osl-3.0.php  * http://www.opensource.org/licenses/osl-3.0.php
11   *  *
12   *******************************************************/  *******************************************************/
13    
14    
15  /* PCG iterations */  /* PCG iterations */
16    
# Line 78  err_t Paso_Solver_PCG( Line 77  err_t Paso_Solver_PCG(
77      Paso_Performance* pp) {      Paso_Performance* pp) {
78    
79    /* Local variables */    /* Local variables */
80    dim_t num_iter=0,maxit,num_iter_global, chunk_size=-1, len,rest, n_chunks, np, ipp;    dim_t num_iter=0,maxit,num_iter_global, chunk_size=-1, len,rest, np, ipp;
81    register double ss,ss1;    register double ss,ss1;
82    dim_t i0, istart, iend;    dim_t i0, istart, iend;
83    bool_t breakFlag=FALSE, maxIterFlag=FALSE, convergeFlag=FALSE;    bool_t breakFlag=FALSE, maxIterFlag=FALSE, convergeFlag=FALSE;
# Line 86  err_t Paso_Solver_PCG( Line 85  err_t Paso_Solver_PCG(
85    dim_t n = Paso_SystemMatrix_getTotalNumRows(A);    dim_t n = Paso_SystemMatrix_getTotalNumRows(A);
86    double *resid = tolerance, *rs=NULL, *p=NULL, *v=NULL, *x2=NULL ;    double *resid = tolerance, *rs=NULL, *p=NULL, *v=NULL, *x2=NULL ;
87    double tau_old,tau,beta,delta,gamma_1,gamma_2,alpha,sum_1,sum_2,sum_3,sum_4,sum_5,tol;    double tau_old,tau,beta,delta,gamma_1,gamma_2,alpha,sum_1,sum_2,sum_3,sum_4,sum_5,tol;
88    double norm_of_residual,norm_of_residual_global, loc_sum[2], sum[2];  #ifdef PASO_MPI
89    register double r_tmp,d,rs_tmp,x2_tmp,x_tmp;    double loc_sum[2], sum[2];
90    char* chksz_chr;  #endif
91    np=omp_get_max_threads();    double norm_of_residual,norm_of_residual_global;
92      register double d;
93    
94      /* Should not be any executable code before this ifdef */
95    
96  #ifdef USE_DYNAMIC_SCHEDULING  #ifdef USE_DYNAMIC_SCHEDULING
97    
98        /* Watch out for these declarations (see above) */
99        char* chksz_chr;
100        dim_t n_chunks;
101    
102      chksz_chr=getenv("PASO_CHUNK_SIZE_PCG");      chksz_chr=getenv("PASO_CHUNK_SIZE_PCG");
103      if (chksz_chr!=NULL) sscanf(chksz_chr, "%d",&chunk_size);      if (chksz_chr!=NULL) sscanf(chksz_chr, "%d",&chunk_size);
104        np=omp_get_max_threads();
105      chunk_size=MIN(MAX(1,chunk_size),n/np);      chunk_size=MIN(MAX(1,chunk_size),n/np);
106      n_chunks=n/chunk_size;      n_chunks=n/chunk_size;
107      if (n_chunks*chunk_size<n) n_chunks+=1;      if (n_chunks*chunk_size<n) n_chunks+=1;
108  #else  #else
109        np=omp_get_max_threads();
110      len=n/np;      len=n/np;
111      rest=n-len*np;      rest=n-len*np;
112  #endif  #endif
# Line 151  err_t Paso_Solver_PCG( Line 160  err_t Paso_Solver_PCG(
160         #endif         #endif
161      }      }
162      num_iter=0;      num_iter=0;
163    
164        /* PGH */
165        /* without this we get a use of an unititialised var below */
166        tau = 0;
167    
168      /* start of iteration */      /* start of iteration */
169      while (!(convergeFlag || maxIterFlag || breakFlag)) {      while (!(convergeFlag || maxIterFlag || breakFlag)) {
170             ++(num_iter);             ++(num_iter);
171    
172               /* PGH */
173               /* The next lines were commented out before I got here */
174             /* v=prec(r)  */             /* v=prec(r)  */
175               /* tau=v*r; */
176               /* leading to the use of an unititialised var below */
177    
178             Performance_stopMonitor(pp,PERFORMANCE_SOLVER);             Performance_stopMonitor(pp,PERFORMANCE_SOLVER);
179             Performance_startMonitor(pp,PERFORMANCE_PRECONDITIONER);             Performance_startMonitor(pp,PERFORMANCE_PRECONDITIONER);
180             Paso_Solver_solvePreconditioner(A,v,r);             Paso_Solver_solvePreconditioner(A,v,r);
181             Performance_stopMonitor(pp,PERFORMANCE_PRECONDITIONER);             Performance_stopMonitor(pp,PERFORMANCE_PRECONDITIONER);
182             Performance_startMonitor(pp,PERFORMANCE_SOLVER);             Performance_startMonitor(pp,PERFORMANCE_SOLVER);
183             /* tau=v*r    */  
184         sum_1 = 0;         sum_1 = 0;
185             #pragma omp parallel private(i0, istart, iend, ipp, ss)             #pragma omp parallel private(i0, istart, iend, ipp, ss)
186             {             {

Legend:
Removed from v.1579  
changed lines
  Added in v.1811

  ViewVC Help
Powered by ViewVC 1.1.26