[cig-commits] r1405 - in branches/s-wang2: for_deal.II/examples/step-32 for_deal.II/include/deal.II/lac for_deal.II/source/lac include/aspect source source/simulator tests
s-wang at dealii.org
s-wang at dealii.org
Fri Nov 30 21:41:37 PST 2012
Author: s-wang
Date: 2012-11-30 22:41:37 -0700 (Fri, 30 Nov 2012)
New Revision: 1405
Added:
branches/s-wang2/include/aspect/global_petsc.h
Modified:
branches/s-wang2/for_deal.II/examples/step-32/test-step-32.cc
branches/s-wang2/for_deal.II/include/deal.II/lac/trilinos_vector_base.h
branches/s-wang2/for_deal.II/source/lac/petsc_precondition.cc
branches/s-wang2/include/aspect/global.h
branches/s-wang2/include/aspect/global_trilinos.h
branches/s-wang2/source/adiabatic_conditions.cc
branches/s-wang2/source/main.cc
branches/s-wang2/source/simulator/assembly.cc
branches/s-wang2/source/simulator/core.cc
branches/s-wang2/source/simulator/helper_functions.cc
branches/s-wang2/source/simulator/initial_conditions.cc
branches/s-wang2/source/simulator/solver.cc
branches/s-wang2/tests/Makefile
Log:
added option in global.h to switch between petsc & trilinos.
Modified: branches/s-wang2/for_deal.II/examples/step-32/test-step-32.cc
===================================================================
--- branches/s-wang2/for_deal.II/examples/step-32/test-step-32.cc 2012-11-30 23:34:02 UTC (rev 1404)
+++ branches/s-wang2/for_deal.II/examples/step-32/test-step-32.cc 2012-12-01 05:41:37 UTC (rev 1405)
@@ -3193,8 +3193,8 @@
Assembly::CopyData::
StokesSystem<dim> (stokes_fe));
- stokes_matrix.compress(dealii::VectorOperation::add);
- stokes_rhs.compress(dealii::VectorOperation::add); //stokes_rhs.compress(Add);
+ stokes_matrix.compress(LinearAlgebra::Add);
+ stokes_rhs.compress(LinearAlgebra::Add); //stokes_rhs.compress(Add);
// std::ofstream filename_matrix, filename_vector;
// filename_matrix.open("stokes_matrix00.txt");
Modified: branches/s-wang2/for_deal.II/include/deal.II/lac/trilinos_vector_base.h
===================================================================
--- branches/s-wang2/for_deal.II/include/deal.II/lac/trilinos_vector_base.h 2012-11-30 23:34:02 UTC (rev 1404)
+++ branches/s-wang2/for_deal.II/include/deal.II/lac/trilinos_vector_base.h 2012-12-01 05:41:37 UTC (rev 1405)
@@ -905,6 +905,8 @@
* to the dealii::Vector<number>
* class.
*/
+
+ void update_ghost_values() const {}; // shuqiangwang
void print (const char *format = 0) const;
/**
Modified: branches/s-wang2/for_deal.II/source/lac/petsc_precondition.cc
===================================================================
--- branches/s-wang2/for_deal.II/source/lac/petsc_precondition.cc 2012-11-30 23:34:02 UTC (rev 1404)
+++ branches/s-wang2/for_deal.II/source/lac/petsc_precondition.cc 2012-12-01 05:41:37 UTC (rev 1405)
@@ -379,16 +379,15 @@
create_pc();
int ierr;
-// ierr = PCSetType (pc, const_cast<char *>(PCILU));
-// ierr = PCSetType (pc, PCSOR);
-// ierr = PCSetType (pc, PCILU);
-// AssertThrow (ierr == 0, ExcPETScError(ierr));
-
- ierr = PCSetType (pc, const_cast<char *>(PCHYPRE));
+ ierr = PCSetType (pc, const_cast<char *>(PCILU));
+ ierr = PCSetType (pc, PCSOR);
AssertThrow (ierr == 0, ExcPETScError(ierr));
- ierr = PCHYPRESetType(pc, "pilut");
- AssertThrow (ierr == 0, ExcPETScError(ierr));
+// ierr = PCSetType (pc, const_cast<char *>(PCHYPRE));
+// AssertThrow (ierr == 0, ExcPETScError(ierr));
+//
+// ierr = PCHYPRESetType(pc, "pilut");
+// AssertThrow (ierr == 0, ExcPETScError(ierr));
// then set flags
PCFactorSetLevels (pc, additional_data.levels);
Modified: branches/s-wang2/include/aspect/global.h
===================================================================
--- branches/s-wang2/include/aspect/global.h 2012-11-30 23:34:02 UTC (rev 1404)
+++ branches/s-wang2/include/aspect/global.h 2012-12-01 05:41:37 UTC (rev 1405)
@@ -23,107 +23,19 @@
#ifndef __aspect__global_h
#define __aspect__global_h
+#define USE_PETSC_LA
+//#define USE_TRILINOS_LA
-#include <deal.II/lac/petsc_parallel_block_vector.h>
-#include <deal.II/lac/petsc_parallel_block_sparse_matrix.h>
-#include <deal.II/lac/petsc_solver.h>
-#include <deal.II/lac/petsc_precondition.h>
-#include <boost/archive/binary_oarchive.hpp>
-#include <boost/archive/binary_iarchive.hpp>
-#include <boost/archive/text_oarchive.hpp>
-#include <boost/archive/text_iarchive.hpp>
-namespace aspect
-{
- /**
- * A variable whose value denotes the number of seconds in one year.
- */
- extern const double year_in_seconds;
+#ifdef USE_PETSC_LA
+#undef USE_TRILINOS
+#include <aspect/global_petsc.h>
+#endif
- /**
- * A variable that denotes whether we should periodically
- * output statistics about memory consumption, run times, etc
- * via the Simulator::output_statistics() function or other
- * means.
- */
- extern const bool output_parallel_statistics;
+#ifdef USE_TRILINOS_LA
+#undef USE_PETSC_LA
+#include <aspect/global_trilinos.h>
+#endif
- /**
- * A typedef that denotes the BOOST stream type for reading data
- * during serialization. The type chosen here is a binary archive
- * which we subsequently will have to un-compress.
- */
- typedef boost::archive::binary_iarchive iarchive;
-
- /**
- * A typedef that denotes the BOOST stream type for writing data
- * during serialization. The type chosen here is a binary archive
- * which we compress before writing it into a file.
- */
- typedef boost::archive::binary_oarchive oarchive;
-
- /**
- * A namespace that contains typedefs for classes used in
- * the linear algebra description.
- */
- namespace LinearAlgebra
- {
- using namespace dealii;
-
-
- /**
- * Typedef for the vector type used.
- */
- typedef PETScWrappers::MPI::Vector Vector;
-
- /**
- * Typedef for the type used to describe vectors that
- * consist of multiple blocks.
- */
- typedef PETScWrappers::MPI::BlockVector BlockVector;
-
- /**
- * Typedef for the sparse matrix type used.
- */
- typedef PETScWrappers::MPI::SparseMatrix SparseMatrix;
-
- /**
- * Typedef for the type used to describe sparse matrices that
- * consist of multiple blocks.
- */
- typedef PETScWrappers::MPI::BlockSparseMatrix BlockSparseMatrix;
-
-// typedef PETScWrappers::SolverCG SolverCG;
-
- /**
- * Typedef for the AMG preconditioner type used for the
- * top left block of the Stokes matrix.
- */
- typedef PETScWrappers::PreconditionBoomerAMG PreconditionAMG;
-
- /**
- * Typedef for the Incomplete Cholesky preconditioner used
- * for other blocks of the system matrix.
- */
- typedef PETScWrappers::PreconditionICC PreconditionIC;
-
- /**
- * Typedef for the Incomplete LU decomposition preconditioner used
- * for other blocks of the system matrix.
- */
- typedef PETScWrappers::PreconditionILU PreconditionILU;
- }
-}
-
-/**
- * A macro that is used in instantiating the ASPECT classes and functions
- * for both 2d and 3d. Call this macro with the name of another macro that
- * when called with a single integer argument instantiates the respective
- * classes in the given space dimension.
- */
-#define ASPECT_INSTANTIATE(INSTANTIATIONS) \
- INSTANTIATIONS(2) \
- INSTANTIATIONS(3)
-
#endif
Added: branches/s-wang2/include/aspect/global_petsc.h
===================================================================
--- branches/s-wang2/include/aspect/global_petsc.h (rev 0)
+++ branches/s-wang2/include/aspect/global_petsc.h 2012-12-01 05:41:37 UTC (rev 1405)
@@ -0,0 +1,135 @@
+/*
+ Copyright (C) 2011, 2012 by the authors of the ASPECT code.
+
+ This file is part of ASPECT.
+
+ ASPECT is free software; you can redistribute it and/or modify
+ it under the terms of the GNU General Public License as published by
+ the Free Software Foundation; either version 2, or (at your option)
+ any later version.
+
+ ASPECT is distributed in the hope that it will be useful,
+ but WITHOUT ANY WARRANTY; without even the implied warranty of
+ MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
+ GNU General Public License for more details.
+
+ You should have received a copy of the GNU General Public License
+ along with ASPECT; see the file doc/COPYING. If not see
+ <http://www.gnu.org/licenses/>.
+*/
+/* $Id: global.h 1394 2012-11-29 05:55:56Z s-wang $ */
+
+
+#ifndef __aspect__global_petsc_h
+#define __aspect__global_petsc_h
+
+
+#include <deal.II/lac/petsc_parallel_block_vector.h>
+#include <deal.II/lac/petsc_parallel_block_sparse_matrix.h>
+#include <deal.II/lac/petsc_precondition.h>
+#include <deal.II/lac/petsc_solver.h>
+
+#include <boost/archive/binary_oarchive.hpp>
+#include <boost/archive/binary_iarchive.hpp>
+#include <boost/archive/text_oarchive.hpp>
+#include <boost/archive/text_iarchive.hpp>
+namespace aspect
+{
+ /**
+ * A variable whose value denotes the number of seconds in one year.
+ */
+ extern const double year_in_seconds;
+
+ /**
+ * A variable that denotes whether we should periodically
+ * output statistics about memory consumption, run times, etc
+ * via the Simulator::output_statistics() function or other
+ * means.
+ */
+ extern const bool output_parallel_statistics;
+
+
+ /**
+ * A typedef that denotes the BOOST stream type for reading data
+ * during serialization. The type chosen here is a binary archive
+ * which we subsequently will have to un-compress.
+ */
+ typedef boost::archive::binary_iarchive iarchive;
+
+ /**
+ * A typedef that denotes the BOOST stream type for writing data
+ * during serialization. The type chosen here is a binary archive
+ * which we compress before writing it into a file.
+ */
+ typedef boost::archive::binary_oarchive oarchive;
+
+ /**
+ * A namespace that contains typedefs for classes used in
+ * the linear algebra description.
+ */
+ namespace LinearAlgebra
+ {
+ using namespace dealii;
+
+
+ /**
+ * Typedef for the vector type used.
+ */
+ typedef PETScWrappers::MPI::Vector Vector;
+
+ /**
+ * Typedef for the type used to describe vectors that
+ * consist of multiple blocks.
+ */
+ typedef PETScWrappers::MPI::BlockVector BlockVector;
+
+ /**
+ * Typedef for the sparse matrix type used.
+ */
+ typedef PETScWrappers::MPI::SparseMatrix SparseMatrix;
+
+ /**
+ * Typedef for the type used to describe sparse matrices that
+ * consist of multiple blocks.
+ */
+ typedef PETScWrappers::MPI::BlockSparseMatrix BlockSparseMatrix;
+
+// typedef PETScWrappers::SolverCG SolverCG;
+
+ /**
+ * Typedef for the AMG preconditioner type used for the
+ * top left block of the Stokes matrix.
+ */
+ typedef PETScWrappers::PreconditionBoomerAMG PreconditionAMG;
+
+ /**
+ * Typedef for the Incomplete Cholesky preconditioner used
+ * for other blocks of the system matrix.
+ */
+ typedef PETScWrappers::PreconditionICC PreconditionIC;
+
+ /**
+ * Typedef for the Incomplete LU decomposition preconditioner used
+ * for other blocks of the system matrix.
+ */
+ typedef PETScWrappers::PreconditionILU PreconditionILU;
+
+ typedef PETScWrappers::SolverCG SolverCG;
+
+
+ const dealii::VectorOperation::values Add = VectorOperation::add;
+ const dealii::VectorOperation::values Insert = VectorOperation::insert;
+ }
+}
+
+/**
+ * A macro that is used in instantiating the ASPECT classes and functions
+ * for both 2d and 3d. Call this macro with the name of another macro that
+ * when called with a single integer argument instantiates the respective
+ * classes in the given space dimension.
+ */
+#define ASPECT_INSTANTIATE(INSTANTIATIONS) \
+ INSTANTIATIONS(2) \
+ INSTANTIATIONS(3)
+
+#endif
Modified: branches/s-wang2/include/aspect/global_trilinos.h
===================================================================
--- branches/s-wang2/include/aspect/global_trilinos.h 2012-11-30 23:34:02 UTC (rev 1404)
+++ branches/s-wang2/include/aspect/global_trilinos.h 2012-12-01 05:41:37 UTC (rev 1405)
@@ -20,13 +20,14 @@
/* $Id: global.h 895 2012-04-10 12:53:27Z bangerth $ */
-#ifndef __aspect__global_h
-#define __aspect__global_h
+#ifndef __aspect__global_trilinos_h
+#define __aspect__global_trilinos_h
#include <deal.II/lac/trilinos_block_vector.h>
#include <deal.II/lac/trilinos_block_sparse_matrix.h>
#include <deal.II/lac/trilinos_precondition.h>
+#include <deal.II/lac/trilinos_solver.h>
#include <boost/archive/binary_oarchive.hpp>
#include <boost/archive/binary_iarchive.hpp>
@@ -110,6 +111,11 @@
* for other blocks of the system matrix.
*/
typedef TrilinosWrappers::PreconditionILU PreconditionILU;
+
+ typedef TrilinosWrappers::SolverCG SolverCG;
+
+ const Epetra_CombineMode Add = ::Add;
+ const Epetra_CombineMode Insert = ::Insert;
}
}
Modified: branches/s-wang2/source/adiabatic_conditions.cc
===================================================================
--- branches/s-wang2/source/adiabatic_conditions.cc 2012-11-30 23:34:02 UTC (rev 1404)
+++ branches/s-wang2/source/adiabatic_conditions.cc 2012-12-01 05:41:37 UTC (rev 1405)
@@ -75,7 +75,6 @@
// approximation here.
const double density = material_model.density(temperatures[i-1], pressures[i-1],
initial_composition,representative_point);
-// std::cout << density << std::endl;
const double alpha = material_model.thermal_expansion_coefficient(temperatures[i-1], pressures[i-1],
initial_composition,representative_point);
const double cp = material_model.specific_heat(temperatures[i-1], pressures[i-1],
Modified: branches/s-wang2/source/main.cc
===================================================================
--- branches/s-wang2/source/main.cc 2012-11-30 23:34:02 UTC (rev 1404)
+++ branches/s-wang2/source/main.cc 2012-12-01 05:41:37 UTC (rev 1405)
@@ -31,7 +31,6 @@
{
using namespace dealii;
Utilities::MPI::MPI_InitFinalize mpi_initialization(argc, argv);
-// PetscInitialize(&argc,&argv,0,0);
try
{
@@ -173,7 +172,6 @@
dealii::GrowingVectorMemory<dealii::PETScWrappers::MPI::Vector>::release_unused_memory ();
dealii::GrowingVectorMemory<dealii::PETScWrappers::Vector>::release_unused_memory ();
-// PetscFinalize();
return 0;
}
Modified: branches/s-wang2/source/simulator/assembly.cc
===================================================================
--- branches/s-wang2/source/simulator/assembly.cc 2012-11-30 23:34:02 UTC (rev 1404)
+++ branches/s-wang2/source/simulator/assembly.cc 2012-12-01 05:41:37 UTC (rev 1405)
@@ -827,18 +827,20 @@
Mp_preconditioner.reset (new LinearAlgebra::PreconditionILU());
Amg_preconditioner.reset (new LinearAlgebra::PreconditionAMG());
-
+#ifdef USE_PETSC_LA
LinearAlgebra::PreconditionAMG::AdditionalData Amg_data(true);
- //Amg_data.constant_modes = constant_modes;
- //Amg_data.elliptic = true;
- //Amg_data.higher_order_elements = true;
- //Amg_data.smoother_sweeps = 2;
- //Amg_data.aggregation_threshold = 0.02;
+#else
+ LinearAlgebra::PreconditionAMG::AdditionalData Amg_data;
+ Amg_data.constant_modes = constant_modes;
+ Amg_data.elliptic = true;
+ Amg_data.higher_order_elements = true;
+ Amg_data.smoother_sweeps = 2;
+ Amg_data.aggregation_threshold = 0.02;
+#endif
Mp_preconditioner->initialize (system_preconditioner_matrix.block(1,1));
Amg_preconditioner->initialize (system_preconditioner_matrix.block(0,0),
Amg_data);
-// exit(0);
rebuild_stokes_preconditioner = false;
@@ -1025,10 +1027,10 @@
StokesSystem<dim> (finite_element));
system_matrix.compress();
- system_rhs.compress(dealii::VectorOperation::add);
+ system_rhs.compress(LinearAlgebra::Add);
if (material_model->is_compressible())
- pressure_shape_function_integrals.compress(dealii::VectorOperation::add);
+ pressure_shape_function_integrals.compress(LinearAlgebra::Add);
rebuild_stokes_matrix = false;
@@ -1047,7 +1049,9 @@
{
preconditioner.reset (new LinearAlgebra::PreconditionILU());
// system_matrix.block(2+index,2+index).write_ascii();
+#ifdef USE_PETSC_LA
if(system_matrix.block(2+index,2+index).linfty_norm()>1e-50) // non-zero
+#endif
preconditioner->initialize (system_matrix.block(2+index,2+index));
}
computing_timer.exit_section();
@@ -1396,7 +1400,7 @@
AdvectionSystem<dim> (finite_element));
system_matrix.compress();
- system_rhs.compress(dealii::VectorOperation::add);
+ system_rhs.compress(LinearAlgebra::Add);
computing_timer.exit_section();
}
Modified: branches/s-wang2/source/simulator/core.cc
===================================================================
--- branches/s-wang2/source/simulator/core.cc 2012-11-30 23:34:02 UTC (rev 1404)
+++ branches/s-wang2/source/simulator/core.cc 2012-12-01 05:41:37 UTC (rev 1405)
@@ -55,6 +55,7 @@
using namespace dealii;
+#ifdef USE_PETSC_LA
/**
* Temporary utility for replacing TrilinosWrappers with PETScWrappers.
*/
@@ -120,8 +121,8 @@
vector.collect_sizes();
}
}
+#endif
-
namespace aspect
{
namespace
@@ -152,7 +153,6 @@
:
parameters (prm),
mpi_communicator (Utilities::MPI::duplicate_communicator (mpi_communicator_)),
-// mpi_communicator(mpi_communicator_),
pcout (std::cout,
(Utilities::MPI::
this_mpi_process(mpi_communicator)
@@ -544,6 +544,7 @@
Simulator<dim>::
setup_system_matrix (const std::vector<IndexSet> &system_partitioning)
{
+ #ifdef USE_TRILINOS_LA
system_matrix.clear ();
TrilinosWrappers::BlockSparsityPattern sp (system_partitioning,
@@ -574,7 +575,8 @@
this_mpi_process(mpi_communicator));
sp.compress();
- //shuqiangwang: this function is not used. system_matrix.reinit (sp);
+ system_matrix.reinit (sp);
+ #endif
}
@@ -583,6 +585,7 @@
void Simulator<dim>::
setup_system_preconditioner (const std::vector<IndexSet> &system_partitioning)
{
+#ifdef USE_TRILINOS_LA
Amg_preconditioner.reset ();
Mp_preconditioner.reset ();
T_preconditioner.reset ();
@@ -608,7 +611,8 @@
this_mpi_process(mpi_communicator));
sp.compress();
- //shuqiangwang; this function is not used. system_preconditioner_matrix.reinit (sp);
+ system_preconditioner_matrix.reinit (sp);
+#endif
}
@@ -783,6 +787,7 @@
constraints.close();
}
+#ifdef USE_PETSC_LA
// finally initialize vectors, matrices, etc.
std::vector<unsigned int> block_sizes, local_sizes;
CIG::convert_block_partitioning(system_dofs_per_block,system_partitioning,block_sizes,local_sizes);
@@ -799,7 +804,23 @@
if (material_model->is_compressible())
pressure_shape_function_integrals.reinit(block_sizes,mpi_communicator,local_sizes); //pressure_shape_function_integrals.reinit (system_partitioning, mpi_communicator);
+#else
+ // finally initialize vectors, matrices, etc.
+ setup_system_matrix (system_partitioning);
+ setup_system_preconditioner (system_partitioning);
+
+ system_rhs.reinit(system_partitioning, mpi_communicator);
+ solution.reinit(system_relevant_partitioning, mpi_communicator);
+ old_solution.reinit(system_relevant_partitioning, mpi_communicator);
+ old_old_solution.reinit(system_relevant_partitioning, mpi_communicator);
+
+ current_linearization_point.reinit (system_relevant_partitioning, MPI_COMM_WORLD);
+
+ if (material_model->is_compressible())
+ pressure_shape_function_integrals.reinit (system_partitioning, mpi_communicator);
+#endif
+
rebuild_stokes_matrix = true;
rebuild_stokes_preconditioner = true;
Modified: branches/s-wang2/source/simulator/helper_functions.cc
===================================================================
--- branches/s-wang2/source/simulator/helper_functions.cc 2012-11-30 23:34:02 UTC (rev 1404)
+++ branches/s-wang2/source/simulator/helper_functions.cc 2012-12-01 05:41:37 UTC (rev 1405)
@@ -78,12 +78,16 @@
<< "* Matrix " << system_matrix.memory_consumption()/mb << std::endl
<< "* 5 Vectors " << 5*solution.memory_consumption()/mb << std::endl
<< "* preconditioner " << (system_preconditioner_matrix.memory_consumption()
- //+ Amg_preconditioner->memory_consumption()
+#ifdef USE_TRILINOS_LA
+ + Amg_preconditioner->memory_consumption()
+#endif
/*+Mp_preconditioner->memory_consumption()
+T_preconditioner->memory_consumption()*/)/mb
<< std::endl
<< " - matrix " << system_preconditioner_matrix.memory_consumption()/mb << std::endl
- // << " - prec vel " << Amg_preconditioner->memory_consumption()/mb << std::endl
+#ifdef USE_TRILINOS_LA
+ << " - prec vel " << Amg_preconditioner->memory_consumption()/mb << std::endl
+#endif
<< " - prec mass " << 0/*Mp_preconditioner->memory_consumption()/mb*/ << std::endl
<< " - prec T " << 0/*T_preconditioner->memory_consumption()/mb*/ << std::endl
<< std::endl;
Modified: branches/s-wang2/source/simulator/initial_conditions.cc
===================================================================
--- branches/s-wang2/source/simulator/initial_conditions.cc 2012-11-30 23:34:02 UTC (rev 1404)
+++ branches/s-wang2/source/simulator/initial_conditions.cc 2012-12-01 05:41:37 UTC (rev 1405)
@@ -129,7 +129,7 @@
}
}
- initial_solution.compress(VectorOperation::insert);
+ initial_solution.compress(LinearAlgebra::Insert);
// we should not have written at all into any of the blocks with
// the exception of the current temperature or composition block
@@ -153,14 +153,8 @@
// then apply constraints and copy the
// result into vectors with ghost elements
-// constraints.print(std::cout);
constraints.distribute(initial_solution);
-// static int debug_index = 0;
-// debug_index++;
-// if(debug_index==2)
-// exit(0);
-
// copy temperature/composition block only
solution.block(2+n) = initial_solution.block(2+n); solution.block(2+n).update_ghost_values();
old_solution.block(2+n) = initial_solution.block(2+n); old_solution.block(2+n).update_ghost_values();
@@ -203,7 +197,6 @@
system_tmp);
system_tmp.compress(); // shuqiangwang: do I need this?
-// system_tmp.print(std::cout,7,false,false);
// we may have hanging nodes, so apply constraints
constraints.distribute (system_tmp);
@@ -294,12 +287,9 @@
cell->set_dof_values (local_projection, system_tmp);
}
system_tmp.compress();
-// system_tmp.print(std::cout,7,true,false);
old_solution.block(1) = system_tmp.block(1);
}
-// old_solution.block(1).print(std::cout,7,true,false);
- //old_solution.compress();
// normalize the pressure in such a way that the surface pressure
// equals a known and desired value
old_solution.update_ghost_values();
@@ -307,8 +297,7 @@
old_solution.update_ghost_values();
// set the current solution to the same value as the previous solution
- solution = old_solution;
- solution.update_ghost_values();
+ solution = old_solution; solution.update_ghost_values();
}
}
Modified: branches/s-wang2/source/simulator/solver.cc
===================================================================
--- branches/s-wang2/source/simulator/solver.cc 2012-11-30 23:34:02 UTC (rev 1404)
+++ branches/s-wang2/source/simulator/solver.cc 2012-12-01 05:41:37 UTC (rev 1405)
@@ -221,9 +221,12 @@
{
SolverControl solver_control(5000, 1e-6 * src.block(1).l2_norm());
+#ifdef USE_PETSC_LA
+ LinearAlgebra::SolverCG solver(solver_control, src.get_mpi_communicator());
+#else
+ LinearAlgebra::SolverCG solver(solver_control);
+#endif
- PETScWrappers::SolverCG solver(solver_control, src.get_mpi_communicator());
-
// Trilinos reports a breakdown
// in case src=dst=0, even
// though it should return
@@ -247,7 +250,11 @@
if (do_solve_A == true)
{
SolverControl solver_control(5000, utmp.l2_norm()*1e-2);
- PETScWrappers::SolverCG solver(solver_control, src.get_mpi_communicator());
+#ifdef USE_PETSC_LA
+ LinearAlgebra::SolverCG solver(solver_control, src.get_mpi_communicator());
+#else
+ LinearAlgebra::SolverCG solver(solver_control);
+#endif
solver.solve(stokes_matrix.block(0,0), dst.block(0), utmp,
a_preconditioner);
}
@@ -286,13 +293,13 @@
// overwrite the vector in residual(), then call set_zero again, and then throw away
// the result
LinearAlgebra::BlockVector
- distributed_solution (system_rhs); //distributed_solution.compress();
+ distributed_solution (system_rhs);
current_constraints.set_zero(distributed_solution); distributed_solution.compress();
// create vector with distribution of system_rhs.
LinearAlgebra::Vector block_remap (system_rhs.block (index+2));
// copy block of current_linearization_point into it, because
// current_linearization is distributed differently.
- block_remap = current_linearization_point.block (index+2); //block_remap.compress();
+ block_remap = current_linearization_point.block (index+2);
// (ab)use the distributed solution vector to temporarily put a residual in
initial_residual = system_matrix.block(index+2,index+2).residual (distributed_solution.block(index+2),
block_remap,
@@ -300,14 +307,18 @@
current_constraints.set_zero(distributed_solution); distributed_solution.compress();
// then overwrite it again with the current best guess and solve the linear system
- distributed_solution.block(index+2) = block_remap; //distributed_solution.compress();
+ distributed_solution.block(index+2) = block_remap;
+#ifdef USE_PETSC_LA
if(system_rhs.block(index+2).linfty_norm()<1e-50)
{
- distributed_solution.block(index+2) = system_rhs.block(index+2);
if(system_matrix.block(index+2,index+2).linfty_norm()>1e-50)
- pcout << "Simulator<dim>::solve_advection(): something wrong" << std::endl;
+ solver.solve (system_matrix.block(index+2,index+2), distributed_solution.block(index+2),
+ system_rhs.block(index+2), index==0?*T_preconditioner:*C_preconditioner);
+ else
+ distributed_solution.block(index+2) = system_rhs.block(index+2);
}
else
+#endif
solver.solve (system_matrix.block(index+2,index+2), distributed_solution.block(index+2),
system_rhs.block(index+2), index==0?*T_preconditioner:*C_preconditioner);
@@ -375,13 +386,13 @@
// then overwrite it again with the current best guess and solve the linear system
distributed_stokes_solution.block(0) = remap.block(0);
- distributed_stokes_solution.block(1) = remap.block(1); //distributed_stokes_solution.compress();
+ distributed_stokes_solution.block(1) = remap.block(1);
// extract Stokes parts of rhs vector
LinearAlgebra::BlockVector distributed_stokes_rhs;
distributed_stokes_rhs.reinit(system_rhs);
distributed_stokes_rhs.block(0) = system_rhs.block(0);
- distributed_stokes_rhs.block(1) = system_rhs.block(1); //distributed_stokes_rhs.compress();
+ distributed_stokes_rhs.block(1) = system_rhs.block(1);
PrimitiveVectorMemory< LinearAlgebra::BlockVector > mem;
Modified: branches/s-wang2/tests/Makefile
===================================================================
--- branches/s-wang2/tests/Makefile 2012-11-30 23:34:02 UTC (rev 1404)
+++ branches/s-wang2/tests/Makefile 2012-12-01 05:41:37 UTC (rev 1405)
@@ -19,6 +19,7 @@
# run times are reported on lines starting with a '|', so
# simply grep them out before running the diff command
$(tests):
+target:
@if test -d $@ ; then \
for i in `cd $@; echo *` ; do \
cat $@/$$i | egrep -v '^\|' > $@/$$i.notime ; \
More information about the CIG-COMMITS
mailing list