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

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

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

revision 4285 by caltinay, Thu Mar 7 01:08:43 2013 UTC revision 4286 by caltinay, Thu Mar 7 04:28:11 2013 UTC
# Line 185  void Paso_FCT_Solver_initialize(const do Line 185  void Paso_FCT_Solver_initialize(const do
185      Performance_stopMonitor(pp,PERFORMANCE_PRECONDITIONER_INIT);      Performance_stopMonitor(pp,PERFORMANCE_PRECONDITIONER_INIT);
186  }  }
187    
188  /* entry point for update proceedures */  /* entry point for update procedures */
189  err_t Paso_FCT_Solver_update(Paso_FCT_Solver *fct_solver, double* u, double *u_old,  Paso_Options* options, Paso_Performance *pp)  err_t Paso_FCT_Solver_update(Paso_FCT_Solver *fct_solver, double* u, double *u_old,  Paso_Options* options, Paso_Performance *pp)
190  {      {    
191      const index_t method=fct_solver->method;      const index_t method=fct_solver->method;
# Line 268  err_t Paso_FCT_Solver_update_LCN(Paso_FC Line 268  err_t Paso_FCT_Solver_update_LCN(Paso_FC
268        if (options->verbose) printf("Paso_FCT_Solver_update_LCN: convergence after %d Gauss-Seidel steps.\n",sweep_max);        if (options->verbose) printf("Paso_FCT_Solver_update_LCN: convergence after %d Gauss-Seidel steps.\n",sweep_max);
269        errorCode=SOLVER_NO_ERROR;        errorCode=SOLVER_NO_ERROR;
270     } else {     } else {
271        if (options->verbose) printf("Paso_FCT_Solver_update_LCN: Gauss-Seidel failed within %d stesp (rel. tolerance %e).\n",sweep_max,RTOL);        if (options->verbose) printf("Paso_FCT_Solver_update_LCN: Gauss-Seidel failed within %d steps (rel. tolerance %e).\n",sweep_max,RTOL);
272        errorCode= SOLVER_MAXITER_REACHED;        errorCode= SOLVER_MAXITER_REACHED;
273     }     }
274     return errorCode;     return errorCode;
# Line 322  err_t Paso_FCT_Solver_updateNL(Paso_FCT_ Line 322  err_t Paso_FCT_Solver_updateNL(Paso_FCT_
322      Paso_FCT_FluxLimiter_setU_tilda(flux_limiter, b); /* u_tilda = m^{-1} b */      Paso_FCT_FluxLimiter_setU_tilda(flux_limiter, b); /* u_tilda = m^{-1} b */
323      /* u_tilda_connector is completed */      /* u_tilda_connector is completed */
324      /********************************************************************************************************************************************/        /********************************************************************************************************************************************/  
325      /* calculate stopping criterium */      /* calculate stopping criterion */
326      norm_u_tilde=Paso_lsup(n, flux_limiter->u_tilde, flux_limiter->mpi_info);      norm_u_tilde=Paso_lsup(n, flux_limiter->u_tilde, flux_limiter->mpi_info);
327      ATOL= rtol * norm_u_tilde + atol ;      ATOL= rtol * norm_u_tilde + atol ;
328      if (options->verbose) printf("Paso_FCT_Solver_updateNL: iteration starts u_tilda lsup = %e (abs. tol = %e)\n",norm_u_tilde,ATOL);      if (options->verbose) printf("Paso_FCT_Solver_updateNL: iteration starts u_tilda lsup = %e (abs. tol = %e)\n",norm_u_tilde,ATOL);
# Line 375  err_t Paso_FCT_Solver_updateNL(Paso_FCT_ Line 375  err_t Paso_FCT_Solver_updateNL(Paso_FCT_
375            } else {            } else {
376            tol *= MIN(MAX(rate*rate, 1e-2), 0.5);            tol *= MIN(MAX(rate*rate, 1e-2), 0.5);
377            }            }
378            /* use BiCGSTab with jacobi preconditioner ( m - omega * L ) */            /* use BiCGStab with jacobi preconditioner ( m - omega * L ) */
379                Paso_zeroes(n,du);                Paso_zeroes(n,du);
380            errorCode = Paso_Solver_BiCGStab(fctp->iteration_matrix, z, du, &cntIter, &tol, pp);            errorCode = Paso_Solver_BiCGStab(fctp->iteration_matrix, z, du, &cntIter, &tol, pp);
381    
# Line 437  err_t Paso_FCT_Solver_updateNL(Paso_FCT_ Line 437  err_t Paso_FCT_Solver_updateNL(Paso_FCT_
437   *        f_{ij} = (m_{ij} - dt (1-theta) d_{ij}) (u_old[j]-u_old[i]) - (m_{ij} + dt theta d_{ij}) (u[j]-u[i])   *        f_{ij} = (m_{ij} - dt (1-theta) d_{ij}) (u_old[j]-u_old[i]) - (m_{ij} + dt theta d_{ij}) (u[j]-u[i])
438   *           *        
439   *     m=fc->mass matrix   *     m=fc->mass matrix
440   *     d=artifical diffusion matrix = L - K = - fc->iteration matrix - fc->transport matrix (away from main diagonal)   *     d=artificial diffusion matrix = L - K = - fc->iteration matrix - fc->transport matrix (away from main diagonal)
441   *   *
442   *   for CN : theta =0.5   *   for CN : theta =0.5
443   *   for BE : theta = 1.   *   for BE : theta = 1.
# Line 533  void Paso_FCT_setAntiDiffusionFlux_BE(Pa Line 533  void Paso_FCT_setAntiDiffusionFlux_BE(Pa
533    }    }
534  }  }
535    
536  /* special version of the ant-diffusive fluxes for the linear Crank-Nicolson scheme  /* special version of the anti-diffusive fluxes for the linear Crank-Nicolson scheme
537   * in fact this is evaluated for u = 2*u_tilde - u_old which is the predictor   * in fact this is evaluated for u = 2*u_tilde - u_old which is the predictor
538   * of the solution of the stabilized problem at time dt using the forward Euler scheme   * of the solution of the stabilized problem at time dt using the forward Euler scheme
539   *   *

Legend:
Removed from v.4285  
changed lines
  Added in v.4286

  ViewVC Help
Powered by ViewVC 1.1.26