diff --git a/src/common/bart_tools.cc b/src/common/bart_tools.cc
index cd94551fd..8e633fdac 100644
--- a/src/common/bart_tools.cc
+++ b/src/common/bart_tools.cc
@@ -7,12 +7,6 @@
 #include "../equations/even_parity.h"
 #include "../aqdata/aq_lsgc"
 
-/** \brief Function to build finite element for general dimensions specified by 
- * user.
- *
- * \parameter prm A reference to processed ParameterHandler instance
- * \return A raw pointer of finite elements derived from FE_Poly
- */
 template <int dim>
 FE_Poly<TensorProductPolynomials<dim>,dim,dim>* build_finite_element
 (ParameterHandler &prm)
@@ -26,11 +20,6 @@ FE_Poly<TensorProductPolynomials<dim>,dim,dim>* build_finite_element
     return new FE_Q<dim> (p_order);
 }
 
-/** \brief Function to build mesh in calculations for general dimensions
- *
- * \parameter prm A reference to processed ParameterHandler instance
- * \return A shared pointer to MeshGenerator<dim> instance
- */
 template <int dim>
 std_cxx11::shared_ptr<MeshGenerator<dim> > build_mesh (ParameterHandler &prm)
 {
@@ -40,11 +29,6 @@ std_cxx11::shared_ptr<MeshGenerator<dim> > build_mesh (ParameterHandler &prm)
   return mesh_class;
 }
 
-/** \brief Function to build pointer to MaterialProperties class.
- *
- * \parameter prm A reference to processed ParameterHandler instance
- * \return A shared pointer to MaterialProperties instance
- */
 std_cxx11::shared_ptr<MaterialProperties> build_material (ParameterHandler &prm)
 {
   std_cxx11::shared_ptr<MaterialProperties> material_class =
@@ -53,28 +37,7 @@ std_cxx11::shared_ptr<MaterialProperties> build_material (ParameterHandler &prm)
 }
 
 template <int dim>
-std_cxx11::shared_ptr<IterationBase<dim> > build_transport_iteration
-(ParameterHandler &prm,
- const std_cxx11::shared_ptr<MeshGenerator<dim> > msh_ptr,
- const std_cxx11::shared_ptr<AQBase<dim> > aqd_ptr)
-{
-  // in future development this builder will be like other two
-  std_cxx11::shared_ptr<IterationBase<dim> > iteration_class =
-  std_cxx11::shared_ptr<IterationBase<dim> > (new IterationBase<dim>(prm,
-                                                                     msh_ptr,
-                                                                     aqd_ptr));
-  return iteration_class;
-}
-
-/** \brief Build specific transport model
- * 
- * \parameter msh_ptr shared_ptr for MeshGenerator<dim> instance
- * \parameter aqd_ptr shared_ptr for AQBase<dim> instance
- * \parameter mat_ptr shared_ptr for MaterialProperties instance
- * \return a shared pointer to transport model
- */
-template <int dim>
-std_cxx11::shared_ptr<TransportBase<dim> > build_transport_model
+std_cxx11::shared_ptr<EquationBase<dim> > build_transport_model
 (ParameterHandler &prm,
  const std_cxx11::shared_ptr<MeshGenerator<dim> > msh_ptr,
  const std_cxx11::shared_ptr<AQBase<dim> > aqd_ptr,
@@ -83,19 +46,14 @@ std_cxx11::shared_ptr<TransportBase<dim> > build_transport_model
   std::string transport_model_name = prm.get ("transport model");
   AssertThrow (transport_model_name!="none",
                ExcMessage("transport model name incorrect or missing"));
-  std_cxx11::shared_ptr<TransportBase<dim> > transport_class;
+  std_cxx11::shared_ptr<EquationBase<dim> > transport_class;
   if (transport_model_name=="ep")
-    transport_class =
-    std_cxx11::shared_ptr<TransportBase<dim> >
+    equation_class =
+    std_cxx11::shared_ptr<EquationBase<dim> >
     (new EvenParity<dim> (prm, msh_ptr, aqd_ptr, mat_ptr));
-  return transport_class;
+  return equation_class;
 }
 
-/** \brief Function to build angular quadrature for general dimensions
- *
- * \parameter prm A processed ParameterHandler instance
- * \return a shared pointer to specific angular quadrature
- */
 template <int dim>
 std_cxx11::shared_ptr<AQBase<dim> >
 build_aq_model (ParameterHandler &prm)
diff --git a/src/common/bart_tools.h b/src/common/bart_tools.h
index a5cd9e756..c1609c6be 100644
--- a/src/common/bart_tools.h
+++ b/src/common/bart_tools.h
@@ -1,45 +1,73 @@
-#ifndef __bart_builder_h__
-#define __bart_builder_h__
+#ifndef __bart_tools_h__
+#define __bart_tools_h__
 
 #include <deal.II/base/parameter_handler.h>
+#include <deal.II/base/conditional_ostream.h>
 #include <deal.II/fe/fe_poly.h>
 
 #include <string>
 #include <map>
 
 #include "../equations/equation_base.h"
+#include "../aqdata/aq_base.h"
+#include "../mesh/mesh_generator.h"
+#include "../common/material_properties.h"
 
 using namespace dealii;
 
+/** \brief Function to build finite element for general dimensions specified by
+ * user.
+ *
+ * \parameter prm A reference to processed ParameterHandler instance
+ * \return A raw pointer of finite elements derived from FE_Poly
+ */
 template <int dim>
 FE_Poly<TensorProductPolynomials<dim>,dim,dim>*
 build_finite_element (ParameterHandler &prm);
 
+/** \brief Function to build mesh in calculations for general dimensions
+ *
+ * \parameter prm A reference to processed ParameterHandler instance
+ * \return A shared pointer to MeshGenerator<dim> instance
+ */
 template <int dim>
 std_cxx11::shared_ptr<MeshGenerator<dim> >
 build_mesh (ParameterHandler &prm);
 
+/** \brief Function to build pointer to MaterialProperties class.
+ *
+ * \parameter prm A reference to processed ParameterHandler instance
+ * \return A shared pointer to MaterialProperties instance
+ */
 std_cxx11::shared_ptr<MaterialProperties> build_material (ParameterHandler &prm);
 
+/** \brief Build specific transport model
+ *
+ * \parameter msh_ptr shared_ptr for MeshGenerator<dim> instance
+ * \parameter aqd_ptr shared_ptr for AQBase<dim> instance
+ * \parameter mat_ptr shared_ptr for MaterialProperties instance
+ * \return a shared pointer to transport model
+ */
 template <int dim>
-std_cxx11::shared_ptr<IterationBase<dim> >
-build_iterative_solver (ParameterHandler &prm,
-                        const std_cxx11::shared_ptr<MeshGenerator<dim> > msh_ptr,
-                        const std_cxx11::shared_ptr<AQBase<dim> > aqd_ptr);
-
-template <int dim>
-std_cxx11::shared_ptr<TransportBase<dim> >
+std_cxx11::shared_ptr<EquationBase<dim> >
 build_transport_model (ParameterHandler &prm,
                        const std_cxx11::shared_ptr<MeshGenerator<dim> > msh_ptr,
                        const std_cxx11::shared_ptr<AQBase<dim> > aqd_ptr,
                        const std_cxx11::shared_ptr<MaterialProperties> mat_ptr);
 
+/** \brief Function to build angular quadrature for general dimensions
+ *
+ * \parameter prm A processed ParameterHandler instance
+ * \return a shared pointer to specific angular quadrature
+ */
 template <int dim>
 std_cxx11::shared_ptr<AQBase<dim> >
 build_aq_model (ParameterHandler &prm);
 
-ConditionalOStream pout(std::cout,
-                        Utilities::MPI::this_mpi_process(MPI_COMM_WORLD)==0);
+/** \brief ConditionalOStream object used to output things on screen with processor 0
+ */
+ConditionalOStream pout (std::cout,
+                         Utilities::MPI::this_mpi_process(MPI_COMM_WORLD)==0);
 
 void radio (std::string str);
 
@@ -57,4 +85,4 @@ void radio (std::string str, bool boolean);
 
 void radio ();
 
-#endif //__bart_builder_h__
+#endif //__bart_tools_h__
diff --git a/src/common/preconditioner_solver.cc b/src/common/preconditioner_solver.cc
index 14d4b5b90..75129fa41 100644
--- a/src/common/preconditioner_solver.cc
+++ b/src/common/preconditioner_solver.cc
@@ -126,110 +126,108 @@ void PreconditionerSolver::initialize_ho_preconditioners
 void PreconditionerSolver::ho_solve
 (std::vector<PETScWrappers::MPI::SparseMatrix*> &ho_syses,
  std::vector<PETScWrappers::MPI::Vector*> &ho_psis,
- std::vector<PETScWrappers::MPI::Vector*> &ho_rhses)
+ std::vector<PETScWrappers::MPI::Vector*> &ho_rhses,
+ unsigned int i/*component number*/)
 {
   AssertThrow (n_total_ho_vars==ho_syses.size(),
                ExcMessage("num of HO system matrices should be equal to total variable number"));
   AssertThrow (n_total_ho_vars==ho_rhses.size(),
                ExcMessage("num of HO system rhs should be equal to total variable number"));
-  for (unsigned int i=0; i<n_total_ho_vars; ++i)
+  if (ho_linear_solver_name=="cg")
   {
-    if (ho_linear_solver_name=="cg")
-    {
-      PETScWrappers::SolverCG
-      solver (*ho_cn[i], MPI_COMM_WORLD);
-      if (ho_preconditioner_name=="amg")
-        solver.solve (*ho_syses[i],
-                      *ho_psis[i],
-                      *ho_rhses[i],
-                      *(pre_ho_amg)[i]);
-      else if (ho_preconditioner_name=="jacobi")
-        solver.solve (*ho_syses[i],
-                      *ho_psis[i],
-                      *ho_rhses[i],
-                      *(pre_ho_jacobi)[i]);
-      else if (ho_preconditioner_name=="bssor")
-        solver.solve (*ho_syses[i],
-                      *ho_psis[i],
-                      *ho_rhses[i],
-                      *(pre_ho_eisenstat)[i]);
-      else if (ho_preconditioner_name=="parasails")
-        solver.solve (*ho_syses[i],
-                      *ho_psis[i],
-                      *ho_rhses[i],
-                      *(pre_ho_parasails)[i]);
-    }
-    else if (ho_linear_solver_name=="bicgstab")
-    {
-      PETScWrappers::SolverBicgstab
-      solver (*ho_cn[i], MPI_COMM_WORLD);
-      if (ho_preconditioner_name=="amg")
-        solver.solve (*ho_syses[i],
-                      *ho_psis[i],
-                      *ho_rhses[i],
-                      *(pre_ho_amg)[i]);
-      else if (ho_preconditioner_name=="jacobi")
-        solver.solve (*ho_syses[i],
-                      *ho_psis[i],
-                      *ho_rhses[i],
-                      *(pre_ho_jacobi)[i]);
-      else if (ho_preconditioner_name=="bssor")
-        solver.solve (*ho_syses[i],
-                      *ho_psis[i],
-                      *ho_rhses[i],
-                      *(pre_ho_eisenstat)[i]);
-      else if (ho_preconditioner_name=="parasails")
-        solver.solve (*ho_syses[i],
-                      *ho_psis[i],
-                      *ho_rhses[i],
-                      *(pre_ho_parasails)[i]);
-    }
-    else if (ho_linear_solver_name=="gmres")
-    {
-      PETScWrappers::SolverGMRES
-      solver (*ho_cn[i], MPI_COMM_WORLD);
-      if (ho_preconditioner_name=="amg")
-        solver.solve (*ho_syses[i],
-                      *ho_psis[i],
-                      *ho_rhses[i],
-                      *(pre_ho_amg)[i]);
-      else if (ho_preconditioner_name=="jacobi")
-        solver.solve (*ho_syses[i],
-                      *ho_psis[i],
-                      *ho_rhses[i],
-                      *(pre_ho_jacobi)[i]);
-      else if (ho_preconditioner_name=="bssor")
-        solver.solve (*ho_syses[i],
-                      *ho_psis[i],
-                      *ho_rhses[i],
-                      *(pre_ho_eisenstat)[i]);
-      else if (ho_preconditioner_name=="parasails")
-        solver.solve (*ho_syses[i],
-                      *ho_psis[i],
-                      *ho_rhses[i],
-                      *(pre_ho_parasails)[i]);
-    }
-    else// if (linear_solver_name=="direct")
+    PETScWrappers::SolverCG
+    solver (*ho_cn[i], MPI_COMM_WORLD);
+    if (ho_preconditioner_name=="amg")
+      solver.solve (*ho_syses[i],
+                    *ho_psis[i],
+                    *ho_rhses[i],
+                    *(pre_ho_amg)[i]);
+    else if (ho_preconditioner_name=="jacobi")
+      solver.solve (*ho_syses[i],
+                    *ho_psis[i],
+                    *ho_rhses[i],
+                    *(pre_ho_jacobi)[i]);
+    else if (ho_preconditioner_name=="bssor")
+      solver.solve (*ho_syses[i],
+                    *ho_psis[i],
+                    *ho_rhses[i],
+                    *(pre_ho_eisenstat)[i]);
+    else if (ho_preconditioner_name=="parasails")
+      solver.solve (*ho_syses[i],
+                    *ho_psis[i],
+                    *ho_rhses[i],
+                    *(pre_ho_parasails)[i]);
+  }
+  else if (ho_linear_solver_name=="bicgstab")
+  {
+    PETScWrappers::SolverBicgstab
+    solver (*ho_cn[i], MPI_COMM_WORLD);
+    if (ho_preconditioner_name=="amg")
+      solver.solve (*ho_syses[i],
+                    *ho_psis[i],
+                    *ho_rhses[i],
+                    *(pre_ho_amg)[i]);
+    else if (ho_preconditioner_name=="jacobi")
+      solver.solve (*ho_syses[i],
+                    *ho_psis[i],
+                    *ho_rhses[i],
+                    *(pre_ho_jacobi)[i]);
+    else if (ho_preconditioner_name=="bssor")
+      solver.solve (*ho_syses[i],
+                    *ho_psis[i],
+                    *ho_rhses[i],
+                    *(pre_ho_eisenstat)[i]);
+    else if (ho_preconditioner_name=="parasails")
+      solver.solve (*ho_syses[i],
+                    *ho_psis[i],
+                    *ho_rhses[i],
+                    *(pre_ho_parasails)[i]);
+  }
+  else if (ho_linear_solver_name=="gmres")
+  {
+    PETScWrappers::SolverGMRES
+    solver (*ho_cn[i], MPI_COMM_WORLD);
+    if (ho_preconditioner_name=="amg")
+      solver.solve (*ho_syses[i],
+                    *ho_psis[i],
+                    *ho_rhses[i],
+                    *(pre_ho_amg)[i]);
+    else if (ho_preconditioner_name=="jacobi")
+      solver.solve (*ho_syses[i],
+                    *ho_psis[i],
+                    *ho_rhses[i],
+                    *(pre_ho_jacobi)[i]);
+    else if (ho_preconditioner_name=="bssor")
+      solver.solve (*ho_syses[i],
+                    *ho_psis[i],
+                    *ho_rhses[i],
+                    *(pre_ho_eisenstat)[i]);
+    else if (ho_preconditioner_name=="parasails")
+      solver.solve (*ho_syses[i],
+                    *ho_psis[i],
+                    *ho_rhses[i],
+                    *(pre_ho_parasails)[i]);
+  }
+  else// if (linear_solver_name=="direct")
+  {
+    if (!ho_direct_init[i])
     {
-      if (!ho_direct_init[i])
-      {
-        ho_direct[i] = std_cxx11::shared_ptr<PETScWrappers::SparseDirectMUMPS>
-        (new PETScWrappers::SparseDirectMUMPS(*ho_cn[i], MPI_COMM_WORLD));
-        if (transport_model_name=="fo" ||
-            (transport_model_name=="ep" && have_reflective_bc))
-          ho_direct[i]->set_symmetric_mode (false);
-        else
-          ho_direct[i]->set_symmetric_mode (true);
-        ho_direct_init[i] = true;
-      }
-      ho_direct[i]->solve (*ho_syses[i],
-                           *ho_psis[i],
-                           *ho_rhses[i]);
+      ho_direct[i] = std_cxx11::shared_ptr<PETScWrappers::SparseDirectMUMPS>
+      (new PETScWrappers::SparseDirectMUMPS(*ho_cn[i], MPI_COMM_WORLD));
+      if (transport_model_name=="fo" ||
+          (transport_model_name=="ep" && have_reflective_bc))
+        ho_direct[i]->set_symmetric_mode (false);
+      else
+        ho_direct[i]->set_symmetric_mode (true);
+      ho_direct_init[i] = true;
     }
-    // the ho_linear_iters are for reporting linear solver status, test purpose only
-    if (ho_linear_solver_name!="direct")
-      ho_linear_iters[i] = ho_cn[i]->last_step ();
+    ho_direct[i]->solve (*ho_syses[i],
+                         *ho_psis[i],
+                         *ho_rhses[i]);
   }
+  // the ho_linear_iters are for reporting linear solver status, test purpose only
+  if (ho_linear_solver_name!="direct")
+    ho_linear_iters[i] = ho_cn[i]->last_step ();
 }
 
 // the following section is for NDA solving/preconditioning
diff --git a/src/transport/equation_base.h b/src/transport/equation_base.h
index 5c0f4f91e..74a01e82e 100644
--- a/src/transport/equation_base.h
+++ b/src/transport/equation_base.h
@@ -19,7 +19,6 @@
 #include <deal.II/base/parameter_handler.h>
 #include <deal.II/base/index_set.h>
 #include <deal.II/base/utilities.h>
-#include <deal.II/base/conditional_ostream.h>
 
 #include <deal.II/distributed/tria.h>
 
diff --git a/src/transport/even_parity.cc b/src/transport/even_parity.cc
index ecc580f84..ca4278e04 100644
--- a/src/transport/even_parity.cc
+++ b/src/transport/even_parity.cc
@@ -1,11 +1,11 @@
-#include "transport_base.h"
 #include "even_parity.h"
 
 template <int dim>
-EvenParity<dim>::EvenParity (ParameterHandler &prm,
-                             const std_cxx11::shared_ptr<MeshGenerator<dim> > msh_ptr,
-                             const std_cxx11::shared_ptr<AQBase<dim> > aqd_ptr,
-                             const std_cxx11::shared_ptr<MaterialProperties> mat_ptr)
+EvenParity<dim>::EvenParity
+(ParameterHandler &prm,
+ const std_cxx11::shared_ptr<MeshGenerator<dim> > msh_ptr,
+ const std_cxx11::shared_ptr<AQBase<dim> > aqd_ptr,
+ const std_cxx11::shared_ptr<MaterialProperties> mat_ptr)
 :
 TransportBase<dim>(prm, msh_ptr, aqd_ptr, mat_ptr)
 {