/[escript]/trunk/paso/src/TFQMR.c
ViewVC logotype

Diff of /trunk/paso/src/TFQMR.c

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

revision 1778 by gross, Tue Sep 9 07:46:02 2008 UTC revision 3259 by jfenwick, Mon Oct 11 01:48:14 2010 UTC
# Line 1  Line 1 
1    
2  /*******************************************************  /*******************************************************
3   *  *
4   *           Copyright 2003-2007 by ACceSS MNRF  * Copyright (c) 2003-2010 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  /* TFQMR iterations */  /* TFQMR iterations */
16    
# Line 22  Line 23 
23  #include <omp.h>  #include <omp.h>
24  #endif  #endif
25    
26  #ifdef PASO_MPI  #ifdef ESYS_MPI
27  #include <mpi.h>  #include <mpi.h>
28  #endif  #endif
29    
# Line 40  Line 41 
41  *  =========  *  =========
42  *  *
43  *  r       (input) DOUBLE PRECISION array, dimension N.  *  r       (input) DOUBLE PRECISION array, dimension N.
44  *          On entry, residual of inital guess x  *        
45  *  *
46  *  x       (input/output) DOUBLE PRECISION array, dimension N.  *  x       (input/output) DOUBLE PRECISION array, dimension N.
47  *          On input, the initial guess.  *        
48  *  *
49  *  ITER    (input/output) INT  *  ITER    (input/output) INT
50  *          On input, the maximum iterations to be performed.  *          On input, the maximum iterations to be performed.
# Line 83  err_t Paso_Solver_TFQMR( Line 84  err_t Paso_Solver_TFQMR(
84    bool_t breakFlag=FALSE, maxIterFlag=FALSE, convergeFlag=FALSE;    bool_t breakFlag=FALSE, maxIterFlag=FALSE, convergeFlag=FALSE;
85    err_t status = SOLVER_NO_ERROR;    err_t status = SOLVER_NO_ERROR;
86    dim_t n = Paso_SystemMatrix_getTotalNumRows(A);    dim_t n = Paso_SystemMatrix_getTotalNumRows(A);
87    double  *u1=NULL, *u2=NULL, *y1=NULL, *y2=NULL, *d=NULL, *w=NULL, *v=NULL, *v_old=NULL;    double  *u1=NULL, *u2=NULL, *y1=NULL, *y2=NULL, *d=NULL, *w=NULL, *v=NULL, *v_old=NULL, *temp_vector=NULL,*res=NULL;
88    
89    double eta,theta,tau,rho,beta,alpha,sigma,rhon,c;    double eta,theta,tau,rho,beta,alpha,sigma,rhon,c;
90    
# Line 104  err_t Paso_Solver_TFQMR( Line 105  err_t Paso_Solver_TFQMR(
105    w=TMPMEMALLOC(n,double);    w=TMPMEMALLOC(n,double);
106    v=TMPMEMALLOC(n,double);    v=TMPMEMALLOC(n,double);
107    v_old=TMPMEMALLOC(n,double);    v_old=TMPMEMALLOC(n,double);
108        temp_vector=TMPMEMALLOC(n,double);
109      res=TMPMEMALLOC(n,double);
110    
111   if (u1 ==NULL || u2== NULL || y1 == NULL || y2== NULL || d==NULL || w==NULL || v==NULL || v_old==NULL) {   if (u1 ==NULL || u2== NULL || y1 == NULL || y2== NULL || d==NULL || w==NULL || v==NULL || v_old==NULL) {
112       status=SOLVER_MEMORY_ERROR;       status=SOLVER_MEMORY_ERROR;
# Line 117  err_t Paso_Solver_TFQMR( Line 119  err_t Paso_Solver_TFQMR(
119      status=SOLVER_INPUT_ERROR;      status=SOLVER_INPUT_ERROR;
120    }    }
121        
122      Paso_zeroes(n,x);
123      
124    Performance_startMonitor(pp,PERFORMANCE_PRECONDITIONER);    Performance_startMonitor(pp,PERFORMANCE_PRECONDITIONER);
125    Paso_Solver_solvePreconditioner(A,r,r);    Paso_SystemMatrix_solvePreconditioner(A,res,r);
126    Performance_stopMonitor(pp,PERFORMANCE_PRECONDITIONER);    Performance_stopMonitor(pp,PERFORMANCE_PRECONDITIONER);
127        
128    Performance_startMonitor(pp,PERFORMANCE_SOLVER);    Performance_startMonitor(pp,PERFORMANCE_SOLVER);
# Line 126  err_t Paso_Solver_TFQMR( Line 130  err_t Paso_Solver_TFQMR(
130    Paso_zeroes(n,u2);    Paso_zeroes(n,u2);
131    Paso_zeroes(n,y2);    Paso_zeroes(n,y2);
132        
133    Paso_Copy(n,w,r);    Paso_Copy(n,w,res);
134    Paso_Copy(n,y1,r);    Paso_Copy(n,y1,res);
135                
136    Paso_zeroes(n,d);    Paso_zeroes(n,d);
137        
138    Performance_stopMonitor(pp,PERFORMANCE_SOLVER);    Performance_stopMonitor(pp,PERFORMANCE_SOLVER);
139    Performance_startMonitor(pp,PERFORMANCE_MVM);    Performance_startMonitor(pp,PERFORMANCE_MVM);
140    Paso_SystemMatrix_MatrixVector_CSR_OFFSET0(ONE, A, y1,ZERO,v);    Paso_SystemMatrix_MatrixVector_CSR_OFFSET0(PASO_ONE, A, y1,PASO_ZERO,temp_vector);
141    Performance_stopMonitor(pp,PERFORMANCE_MVM);    Performance_stopMonitor(pp,PERFORMANCE_MVM);
142    Performance_startMonitor(pp,PERFORMANCE_SOLVER);    Performance_startMonitor(pp,PERFORMANCE_SOLVER);
143        
144    Performance_stopMonitor(pp,PERFORMANCE_SOLVER);    Performance_stopMonitor(pp,PERFORMANCE_SOLVER);
145    Performance_startMonitor(pp,PERFORMANCE_PRECONDITIONER);    Performance_startMonitor(pp,PERFORMANCE_PRECONDITIONER);
146    Paso_Solver_solvePreconditioner(A,v,v);    Paso_SystemMatrix_solvePreconditioner(A,v,temp_vector);
147    Performance_stopMonitor(pp,PERFORMANCE_PRECONDITIONER);    Performance_stopMonitor(pp,PERFORMANCE_PRECONDITIONER);
148    Performance_startMonitor(pp,PERFORMANCE_SOLVER);    Performance_startMonitor(pp,PERFORMANCE_SOLVER);
149        
# Line 148  err_t Paso_Solver_TFQMR( Line 152  err_t Paso_Solver_TFQMR(
152    theta = 0.0;    theta = 0.0;
153    eta = 0.0;    eta = 0.0;
154        
155    tau = Paso_l2(n,r,A->mpi_info);    tau = Paso_l2(n,res,A->mpi_info);
156        
157    rho = tau * tau;    rho = tau * tau;
158                
# Line 158  err_t Paso_Solver_TFQMR( Line 162  err_t Paso_Solver_TFQMR(
162    {    {
163                        
164    
165       sigma=Paso_InnerProduct(n,r,v,A->mpi_info);       sigma=Paso_InnerProduct(n,res,v,A->mpi_info);
166            
167       if (! (breakFlag = (ABS(sigma) == 0.))) {       if (! (breakFlag = (ABS(sigma) == 0.))) {
168            
# Line 172  err_t Paso_Solver_TFQMR( Line 176  err_t Paso_Solver_TFQMR(
176                        
177            Performance_stopMonitor(pp,PERFORMANCE_SOLVER);            Performance_stopMonitor(pp,PERFORMANCE_SOLVER);
178            Performance_startMonitor(pp,PERFORMANCE_MVM);            Performance_startMonitor(pp,PERFORMANCE_MVM);
179            Paso_SystemMatrix_MatrixVector_CSR_OFFSET0(ONE, A, y2,ZERO,u2);            Paso_SystemMatrix_MatrixVector_CSR_OFFSET0(PASO_ONE, A, y2,PASO_ZERO,temp_vector);
180            Performance_stopMonitor(pp,PERFORMANCE_MVM);            Performance_stopMonitor(pp,PERFORMANCE_MVM);
181            Performance_startMonitor(pp,PERFORMANCE_SOLVER);            Performance_startMonitor(pp,PERFORMANCE_SOLVER);
182                        
183            Performance_stopMonitor(pp,PERFORMANCE_SOLVER);            Performance_stopMonitor(pp,PERFORMANCE_SOLVER);
184            Performance_startMonitor(pp,PERFORMANCE_PRECONDITIONER);            Performance_startMonitor(pp,PERFORMANCE_PRECONDITIONER);
185            Paso_Solver_solvePreconditioner(A,u2,u2);        Paso_SystemMatrix_solvePreconditioner(A,u2,temp_vector);
186            Performance_stopMonitor(pp,PERFORMANCE_PRECONDITIONER);            Performance_stopMonitor(pp,PERFORMANCE_PRECONDITIONER);
187            Performance_startMonitor(pp,PERFORMANCE_SOLVER);            Performance_startMonitor(pp,PERFORMANCE_SOLVER);
188           }           }
# Line 201  err_t Paso_Solver_TFQMR( Line 205  err_t Paso_Solver_TFQMR(
205    
206       breakFlag = (ABS(rho) == 0);       breakFlag = (ABS(rho) == 0);
207    
208       rhon = Paso_InnerProduct(n,r,w,A->mpi_info);       rhon = Paso_InnerProduct(n,res,w,A->mpi_info);
209       beta = rhon / rho;       beta = rhon / rho;
210       rho = rhon;       rho = rhon;
211        
# Line 210  err_t Paso_Solver_TFQMR( Line 214  err_t Paso_Solver_TFQMR(
214       Performance_stopMonitor(pp,PERFORMANCE_SOLVER);       Performance_stopMonitor(pp,PERFORMANCE_SOLVER);
215            
216       Performance_startMonitor(pp,PERFORMANCE_MVM);       Performance_startMonitor(pp,PERFORMANCE_MVM);
217       Paso_SystemMatrix_MatrixVector_CSR_OFFSET0(ONE, A, y1,ZERO,u1);       Paso_SystemMatrix_MatrixVector_CSR_OFFSET0(PASO_ONE, A, y1,PASO_ZERO,temp_vector);
218       Performance_stopMonitor(pp,PERFORMANCE_MVM);       Performance_stopMonitor(pp,PERFORMANCE_MVM);
219            
220       Performance_startMonitor(pp,PERFORMANCE_PRECONDITIONER);       Performance_startMonitor(pp,PERFORMANCE_PRECONDITIONER);
221       Paso_Solver_solvePreconditioner(A,u1,u1);       Paso_SystemMatrix_solvePreconditioner(A,u1,temp_vector);
222       Performance_stopMonitor(pp,PERFORMANCE_PRECONDITIONER);       Performance_stopMonitor(pp,PERFORMANCE_PRECONDITIONER);
223            
224       Performance_startMonitor(pp,PERFORMANCE_SOLVER);       Performance_startMonitor(pp,PERFORMANCE_SOLVER);
# Line 247  err_t Paso_Solver_TFQMR( Line 251  err_t Paso_Solver_TFQMR(
251      TMPMEMFREE(w);      TMPMEMFREE(w);
252      TMPMEMFREE(v);      TMPMEMFREE(v);
253      TMPMEMFREE(v_old);      TMPMEMFREE(v_old);
254          TMPMEMFREE(temp_vector);
255        TMPMEMFREE(res);
256      *iter=num_iter;      *iter=num_iter;
257      *tolerance=norm_of_residual;      *tolerance=norm_of_residual;
258            

Legend:
Removed from v.1778  
changed lines
  Added in v.3259

  ViewVC Help
Powered by ViewVC 1.1.26