C++ Adapt::set_error_form方法代码示例

int main(int argc, char* argv[])
{
  // Time measurement.
  TimePeriod cpu_time;
  cpu_time.tick();

  // Load the mesh.
  Mesh u_mesh, v_mesh;
  H2DReader mloader;
  mloader.load("bracket.mesh", &u_mesh);

  // Initial mesh refinements.
  u_mesh.refine_element_id(1);
  u_mesh.refine_element_id(4);

  // Create initial mesh for the vertical displacement component.
  // This also initializes the multimesh hp-FEM.
  v_mesh.copy(&u_mesh);

  // Enter boundary markers.
  BCTypes bc_types;
  bc_types.add_bc_dirichlet(BDY_LEFT);
  bc_types.add_bc_neumann(Hermes::vector<int>(BDY_TOP, BDY_REST));

  // Enter Dirichlet boundary values.
  BCValues bc_values;
  bc_values.add_zero(BDY_LEFT);
    
  // Create H1 spaces with default shapesets.
  H1Space u_space(&u_mesh, &bc_types, &bc_values, P_INIT);
  H1Space v_space(MULTI ? &v_mesh : &u_mesh, &bc_types, &bc_values, P_INIT);

  // Initialize the weak formulation.
  WeakForm wf(2);
  wf.add_matrix_form(0, 0, callback(bilinear_form_0_0), HERMES_SYM);  // note that only one symmetric part is
  wf.add_matrix_form(0, 1, callback(bilinear_form_0_1), HERMES_SYM);  // added in the case of symmetric bilinear
  wf.add_matrix_form(1, 1, callback(bilinear_form_1_1), HERMES_SYM);  // forms
  wf.add_vector_form_surf(1, linear_form_surf_1, linear_form_surf_1_ord, BDY_TOP);

  // Initialize coarse and reference mesh solutions.
  Solution u_sln, v_sln, u_ref_sln, v_ref_sln;

  // Initialize refinement selector.
  H1ProjBasedSelector selector(CAND_LIST, CONV_EXP, H2DRS_DEFAULT_ORDER);

  // DOF and CPU convergence graphs.
  SimpleGraph graph_dof_est, graph_cpu_est;

  // Adaptivity loop:
  int as = 1; 
  bool done = false;
  do
  {
    info("---- Adaptivity step %d:", as);

    // Construct globally refined reference mesh and setup reference space.
    Hermes::vector<Space *>* ref_spaces = construct_refined_spaces(Hermes::vector<Space *>(&u_space, &v_space));

    // Assemble the reference problem.
    info("Solving on reference mesh.");
    bool is_linear = true;
    DiscreteProblem* dp = new DiscreteProblem(&wf, *ref_spaces, is_linear);
    SparseMatrix* matrix = create_matrix(matrix_solver);
    Vector* rhs = create_vector(matrix_solver);
    Solver* solver = create_linear_solver(matrix_solver, matrix, rhs);
    dp->assemble(matrix, rhs);

    // Time measurement.
    cpu_time.tick();
    
    // Solve the linear system of the reference problem. If successful, obtain the solutions.
    if(solver->solve()) Solution::vector_to_solutions(solver->get_solution(), *ref_spaces, 
                                                      Hermes::vector<Solution *>(&u_ref_sln, &v_ref_sln));
    else error ("Matrix solver failed.\n");
  
    // Time measurement.
    cpu_time.tick();

    // Project the fine mesh solution onto the coarse mesh.
    info("Projecting reference solution on coarse mesh.");
    OGProjection::project_global(Hermes::vector<Space *>(&u_space, &v_space), Hermes::vector<Solution *>(&u_ref_sln, &v_ref_sln), 
                   Hermes::vector<Solution *>(&u_sln, &v_sln), matrix_solver); 

    // Calculate element errors.
    info("Calculating error estimate and exact error."); 
    Adapt* adaptivity = new Adapt(Hermes::vector<Space *>(&u_space, &v_space));
    adaptivity->set_error_form(0, 0, bilinear_form_0_0<scalar, scalar>, bilinear_form_0_0<Ord, Ord>);
    adaptivity->set_error_form(0, 1, bilinear_form_0_1<scalar, scalar>, bilinear_form_0_1<Ord, Ord>);
    adaptivity->set_error_form(1, 0, bilinear_form_1_0<scalar, scalar>, bilinear_form_1_0<Ord, Ord>);
    adaptivity->set_error_form(1, 1, bilinear_form_1_1<scalar, scalar>, bilinear_form_1_1<Ord, Ord>);

    // Calculate error estimate for each solution component and the total error estimate.
    Hermes::vector<double> err_est_rel;
    double err_est_rel_total = adaptivity->calc_err_est(Hermes::vector<Solution *>(&u_sln, &v_sln), 
                               Hermes::vector<Solution *>(&u_ref_sln, &v_ref_sln), 
                               &err_est_rel) * 100;

    // Time measurement.
    cpu_time.tick();

//.........这里部分代码省略.........

//.........这里部分代码省略.........
        dp->assemble(matrix, rhs);

        // Time measurement.
        cpu_time.tick();

        // Solve the linear system of the reference problem. If successful, obtain the solutions.
        if(solver->solve()) Solution::vector_to_solutions(solver->get_solution(), *ref_spaces,
                    Tuple<Solution *>(&ref_xdisp_sln, &ref_ydisp_sln, &ref_temp_sln));
        else error ("Matrix solver failed.\n");

        // Time measurement.
        cpu_time.tick();

        // Project the fine mesh solution onto the coarse mesh.
        info("Projecting reference solution on coarse mesh.");
        OGProjection::project_global(Tuple<Space *>(&xdisp, &ydisp, &temp), Tuple<Solution *>(&ref_xdisp_sln, &ref_ydisp_sln, &ref_temp_sln),
                                     Tuple<Solution *>(&xdisp_sln, &ydisp_sln, &temp_sln), matrix_solver);

        // View the coarse mesh solution and polynomial orders.
        s_view_0.show(&xdisp_sln);
        o_view_0.show(&xdisp);
        s_view_1.show(&ydisp_sln);
        o_view_1.show(&ydisp);
        s_view_2.show(&temp_sln);
        o_view_2.show(&temp);

        // Skip visualization time.
        cpu_time.tick(HERMES_SKIP);

        // Calculate element errors.
        info("Calculating error estimate and exact error.");
        Adapt* adaptivity = new Adapt(Tuple<Space *>(&xdisp, &ydisp, &temp),
                                      Tuple<ProjNormType>(HERMES_H1_NORM, HERMES_H1_NORM, HERMES_H1_NORM));
        adaptivity->set_error_form(0, 0, bilinear_form_0_0<scalar, scalar>, bilinear_form_0_0<Ord, Ord>);
        adaptivity->set_error_form(0, 1, bilinear_form_0_1<scalar, scalar>, bilinear_form_0_1<Ord, Ord>);
        adaptivity->set_error_form(0, 2, bilinear_form_0_2<scalar, scalar>, bilinear_form_0_2<Ord, Ord>);
        adaptivity->set_error_form(1, 0, bilinear_form_1_0<scalar, scalar>, bilinear_form_1_0<Ord, Ord>);
        adaptivity->set_error_form(1, 1, bilinear_form_1_1<scalar, scalar>, bilinear_form_1_1<Ord, Ord>);
        adaptivity->set_error_form(1, 2, bilinear_form_1_2<scalar, scalar>, bilinear_form_1_2<Ord, Ord>);
        adaptivity->set_error_form(2, 2, bilinear_form_2_2<scalar, scalar>, bilinear_form_2_2<Ord, Ord>);

        // Calculate error estimate for each solution component and the total error estimate.
        Tuple<double> err_est_rel;
        bool solutions_for_adapt = true;
        double err_est_rel_total = adaptivity->calc_err_est(Tuple<Solution *>(&xdisp_sln, &ydisp_sln, &temp_sln),
                                   Tuple<Solution *>(&ref_xdisp_sln, &ref_ydisp_sln, &ref_temp_sln), solutions_for_adapt,
                                   HERMES_TOTAL_ERROR_REL | HERMES_ELEMENT_ERROR_ABS, &err_est_rel) * 100;

        // Time measurement.
        cpu_time.tick();

        // Report results.
        info("ndof_coarse[xdisp]: %d, ndof_fine[xdisp]: %d, err_est_rel[xdisp]: %g%%",
             xdisp.Space::get_num_dofs(), Space::get_num_dofs((*ref_spaces)[0]), err_est_rel[0]*100);
        info("ndof_coarse[ydisp]: %d, ndof_fine[ydisp]: %d, err_est_rel[ydisp]: %g%%",
             ydisp.Space::get_num_dofs(), Space::get_num_dofs((*ref_spaces)[1]), err_est_rel[1]*100);
        info("ndof_coarse[temp]: %d, ndof_fine[temp]: %d, err_est_rel[temp]: %g%%",
             temp.Space::get_num_dofs(), Space::get_num_dofs((*ref_spaces)[2]), err_est_rel[2]*100);
        info("ndof_coarse_total: %d, ndof_fine_total: %d, err_est_rel_total: %g%%",
             Space::get_num_dofs(Tuple<Space *>(&xdisp, &ydisp, &temp)), Space::get_num_dofs(*ref_spaces), err_est_rel_total);

        // Add entry to DOF and CPU convergence graphs.
        graph_dof_est.add_values(Space::get_num_dofs(Tuple<Space *>(&xdisp, &ydisp, &temp)), err_est_rel_total);
        graph_dof_est.save("conv_dof_est.dat");
        graph_cpu_est.add_values(cpu_time.accumulated(), err_est_rel_total);
        graph_cpu_est.save("conv_cpu_est.dat");

//.........这里部分代码省略.........
    int as = 1; 
    bool done = false;
    do
    {
      info("---- Adaptivity step %d:", as);

      // Construct globally refined reference mesh and setup reference space.
      Tuple<Space *>* ref_spaces = construct_refined_spaces(Tuple<Space *>(&T_space, &M_space));

      // Assemble the reference problem.
      info("Solving on reference mesh.");
      bool is_linear = true;
      DiscreteProblem* dp = new DiscreteProblem(&wf, *ref_spaces, is_linear);
      SparseMatrix* matrix = create_matrix(matrix_solver);
      Vector* rhs = create_vector(matrix_solver);
      Solver* solver = create_linear_solver(matrix_solver, matrix, rhs);
      dp->assemble(matrix, rhs);

      // Now we can deallocate the previous fine meshes.
      if(as > 1){ delete T_fine.get_mesh(); delete M_fine.get_mesh(); }

      // Solve the linear system of the reference problem. If successful, obtain the solutions.
      if(solver->solve()) Solution::vector_to_solutions(solver->get_solution(), *ref_spaces, 
                                              Tuple<Solution *>(&T_fine, &M_fine));
      else error ("Matrix solver failed.\n");

      // Project the fine mesh solution onto the coarse mesh.
      info("Projecting reference solution on coarse mesh.");
      OGProjection::project_global(Tuple<Space *>(&T_space, &M_space), Tuple<Solution *>(&T_fine, &M_fine), 
                     Tuple<Solution *>(&T_coarse, &M_coarse), matrix_solver); 

      // Calculate element errors and total error estimate.
      info("Calculating error estimate."); 
      Adapt* adaptivity = new Adapt(Tuple<Space *>(&T_space, &M_space), Tuple<ProjNormType>(HERMES_H1_NORM, HERMES_H1_NORM));
      adaptivity->set_error_form(0, 0, callback(bilinear_form_sym_0_0));
      adaptivity->set_error_form(0, 1, callback(bilinear_form_sym_0_1));
      adaptivity->set_error_form(1, 0, callback(bilinear_form_sym_1_0));
      adaptivity->set_error_form(1, 1, callback(bilinear_form_sym_1_1));
      double err_est_rel_total = adaptivity->calc_err_est(Tuple<Solution *>(&T_coarse, &M_coarse), 
                                 Tuple<Solution *>(&T_fine, &M_fine), 
                                 HERMES_TOTAL_ERROR_REL | HERMES_ELEMENT_ERROR_ABS) * 100;

      // Report results.
      info("ndof_coarse: %d, ndof_fine: %d, err_est_rel: %g%%", 
        Space::get_num_dofs(Tuple<Space *>(&T_space, &M_space)), Space::get_num_dofs(*ref_spaces), err_est_rel_total);
      
      // Show new coarse meshes and solutions.
      char title[100];
      sprintf(title, "Temperature, t = %g days", CURRENT_TIME/3600./24);
      T_sln_view.set_title(title);
      T_sln_view.show(&T_coarse);
      sprintf(title, "Moisture, t = %g days", CURRENT_TIME/3600./24);
      M_sln_view.set_title(title);
      M_sln_view.show(&M_coarse);
      T_order_view.show(&T_space);
      M_order_view.show(&M_space);

      // If err_est too large, adapt the mesh.
      if (err_est_rel_total < ERR_STOP) 
        done = true;
      else 
      {
        info("Adapting coarse mesh.");
        done = adaptivity->adapt(Tuple<RefinementSelectors::Selector *>(&selector, &selector), 
                                 THRESHOLD, STRATEGY, MESH_REGULARITY);
        if (Space::get_num_dofs(Tuple<Space *>(&T_space, &M_space)) >= NDOF_STOP) 
          done = true;
        else
          // Increase the counter of performed adaptivity steps.
          as++;
      }

      // Clean up.
      delete solver;
      delete matrix;
      delete rhs;
      delete adaptivity;
      delete ref_spaces;
      delete dp;
      
      // Increase counter.
      as++;
    }
    while (done == false);

    // Update time.
    CURRENT_TIME += TAU;

    
    // Save fine mesh solutions for the next time step.
    T_prev.copy(&T_fine);
    M_prev.copy(&M_fine);

    ts++;
  }

  // Wait for all views to be closed.
  View::wait();
  return 0;
}

//.........这里部分代码省略.........
      Hermes::vector<Space *>* ref_spaces = construct_refined_spaces(Hermes::vector<Space *>(&T_space, &M_space));

      // Assemble the reference problem.
      info("Solving on reference mesh.");
      bool is_linear = true;
      DiscreteProblem* dp = new DiscreteProblem(&wf, *ref_spaces, is_linear);
      SparseMatrix* matrix = create_matrix(matrix_solver);
      Vector* rhs = create_vector(matrix_solver);
      Solver* solver = create_linear_solver(matrix_solver, matrix, rhs);
      dp->assemble(matrix, rhs);

      // Now we can deallocate the previous fine meshes.
      if(as > 1){ delete T_fine.get_mesh(); delete M_fine.get_mesh(); }
      
      // Time measurement.
      cpu_time.tick();

      // Solve the linear system of the reference problem. If successful, obtain the solutions.
      if(solver->solve()) Solution::vector_to_solutions(solver->get_solution(), *ref_spaces, 
                                              Hermes::vector<Solution *>(&T_fine, &M_fine));
      else error ("Matrix solver failed.\n");
    
      // Time measurement.
      cpu_time.tick();

      // Project the fine mesh solution onto the coarse mesh.
      info("Projecting reference solution on coarse mesh.");
      OGProjection::project_global(Hermes::vector<Space *>(&T_space, &M_space), Hermes::vector<Solution *>(&T_fine, &M_fine), 
                     Hermes::vector<Solution *>(&T_coarse, &M_coarse), matrix_solver); 

      // Calculate element errors and total error estimate.
      info("Calculating error estimate."); 
      Adapt* adaptivity = new Adapt(Hermes::vector<Space *>(&T_space, &M_space));
      adaptivity->set_error_form(0, 0, callback(bilinear_form_sym_0_0));
      adaptivity->set_error_form(0, 1, callback(bilinear_form_sym_0_1));
      adaptivity->set_error_form(1, 0, callback(bilinear_form_sym_1_0));
      adaptivity->set_error_form(1, 1, callback(bilinear_form_sym_1_1));
      double err_est_rel_total = adaptivity->calc_err_est(Hermes::vector<Solution *>(&T_coarse, &M_coarse), 
                                 Hermes::vector<Solution *>(&T_fine, &M_fine)) * 100;

      // Time measurement.
      cpu_time.tick();

      // Report results.
      info("ndof_coarse: %d, ndof_fine: %d, err_est_rel: %g%%", 
        Space::get_num_dofs(Hermes::vector<Space *>(&T_space, &M_space)), Space::get_num_dofs(*ref_spaces), err_est_rel_total);
      
      // If err_est too large, adapt the mesh.
      if (err_est_rel_total < ERR_STOP) 
        done = true;
      else 
      {
        info("Adapting coarse mesh.");
        done = adaptivity->adapt(Hermes::vector<RefinementSelectors::Selector *>(&selector, &selector), 
                                 THRESHOLD, STRATEGY, MESH_REGULARITY);
        if (Space::get_num_dofs(Hermes::vector<Space *>(&T_space, &M_space)) >= NDOF_STOP) 
          done = true;
        else
          // Increase the counter of performed adaptivity steps.
          as++;
      }

      // Clean up.
      delete solver;
      delete matrix;
      delete rhs;

int main(int argc, char* argv[])
{
  // Time measurement.
  TimePeriod cpu_time;
  cpu_time.tick();

  // Load the mesh.
  Mesh u_mesh, v_mesh;
  H2DReader mloader;
  mloader.load("bracket.mesh", &u_mesh);

  // Initial mesh refinements.
  u_mesh.refine_element(1);
  u_mesh.refine_element(4);

  // Create initial mesh for the vertical displacement component.
  // This also initializes the multimesh hp-FEM.
  v_mesh.copy(&u_mesh);

  // Create H1 spaces with default shapesets.
  H1Space u_space(&u_mesh, bc_types, essential_bc_values, P_INIT);
  H1Space v_space(MULTI ? &v_mesh : &u_mesh, bc_types, essential_bc_values, P_INIT);

  // Initialize the weak formulation.
  WeakForm wf(2);
  wf.add_matrix_form(0, 0, callback(bilinear_form_0_0), HERMES_SYM);  // note that only one symmetric part is
  wf.add_matrix_form(0, 1, callback(bilinear_form_0_1), HERMES_SYM);  // added in the case of symmetric bilinear
  wf.add_matrix_form(1, 1, callback(bilinear_form_1_1), HERMES_SYM);  // forms
  wf.add_vector_form_surf(1, linear_form_surf_1, linear_form_surf_1_ord, BDY_TOP);

  // Initialize coarse and reference mesh solutions.
  Solution u_sln, v_sln, u_ref_sln, v_ref_sln;

  // Initialize refinement selector.
  H1ProjBasedSelector selector(CAND_LIST, CONV_EXP, H2DRS_DEFAULT_ORDER);

  // Initialize views.
  ScalarView s_view_0("Solution[0]", new WinGeom(0, 0, 400, 300));
  s_view_0.show_mesh(false);
  ScalarView s_view_1("Solution[1]", new WinGeom(780, 0, 400, 300));
  s_view_1.show_mesh(false);
  OrderView  o_view_0("Mesh[0]", new WinGeom(410, 0, 360, 300));
  OrderView  o_view_1("Mesh[1]", new WinGeom(1190, 0, 400, 300));

  // DOF and CPU convergence graphs.
  SimpleGraph graph_dof_est, graph_cpu_est;

  // Adaptivity loop:
  int as = 1; 
  bool done = false;
  do
  {
    info("---- Adaptivity step %d:", as);

    // Construct globally refined reference mesh and setup reference space.
    Tuple<Space *>* ref_spaces = construct_refined_spaces(Tuple<Space *>(&u_space, &v_space));

    // Assemble the reference problem.
    info("Solving on reference mesh.");
    bool is_linear = true;
    DiscreteProblem* dp = new DiscreteProblem(&wf, *ref_spaces, is_linear);
    SparseMatrix* matrix = create_matrix(matrix_solver);
    Vector* rhs = create_vector(matrix_solver);
    Solver* solver = create_linear_solver(matrix_solver, matrix, rhs);
    dp->assemble(matrix, rhs);

    // Time measurement.
    cpu_time.tick();
    
    // Solve the linear system of the reference problem. If successful, obtain the solutions.
    if(solver->solve()) Solution::vector_to_solutions(solver->get_solution(), *ref_spaces, 
                                            Tuple<Solution *>(&u_ref_sln, &v_ref_sln));
    else error ("Matrix solver failed.\n");
  
    // Time measurement.
    cpu_time.tick();

    // Project the fine mesh solution onto the coarse mesh.
    info("Projecting reference solution on coarse mesh.");
    OGProjection::project_global(Tuple<Space *>(&u_space, &v_space), Tuple<Solution *>(&u_ref_sln, &v_ref_sln), 
                   Tuple<Solution *>(&u_sln, &v_sln), matrix_solver); 
   
    // View the coarse mesh solution and polynomial orders.
    s_view_0.show(&u_sln); 
    o_view_0.show(&u_space);
    s_view_1.show(&v_sln); 
    o_view_1.show(&v_space);

    // Skip visualization time.
    cpu_time.tick(HERMES_SKIP);

    // Calculate element errors.
    info("Calculating error estimate and exact error."); 
    Adapt* adaptivity = new Adapt(Tuple<Space *>(&u_space, &v_space), Tuple<ProjNormType>(HERMES_H1_NORM, HERMES_H1_NORM));
    adaptivity->set_error_form(0, 0, bilinear_form_0_0<scalar, scalar>, bilinear_form_0_0<Ord, Ord>);
    adaptivity->set_error_form(0, 1, bilinear_form_0_1<scalar, scalar>, bilinear_form_0_1<Ord, Ord>);
    adaptivity->set_error_form(1, 0, bilinear_form_1_0<scalar, scalar>, bilinear_form_1_0<Ord, Ord>);
    adaptivity->set_error_form(1, 1, bilinear_form_1_1<scalar, scalar>, bilinear_form_1_1<Ord, Ord>);

    // Calculate error estimate for each solution component and the total error estimate.
//.........这里部分代码省略.........

//.........这里部分代码省略.........
    DiscreteProblem<std::complex<double> > dp(&wf, ref_space);

    // Time measurement.
    cpu_time.tick();

    // Perform Newton's iteration.
    Hermes::Hermes2D::NewtonSolver<std::complex<double> > newton(&dp);
    try
    {
      newton.set_newton_max_iter(NEWTON_MAX_ITER);
      newton.set_newton_tol(NEWTON_TOL);
      newton.solve();
    }
    catch(Hermes::Exceptions::Exception e)
    {
      e.print_msg();
      throw Hermes::Exceptions::Exception("Newton's iteration failed.");
    };
    // Translate the resulting coefficient vector into the Solution<std::complex<double> > sln.
    Hermes::Hermes2D::Solution<std::complex<double> >::vector_to_solution(newton.get_sln_vector(), ref_space, &ref_sln);
  
    // Project the fine mesh solution onto the coarse mesh.
    Hermes::Mixins::Loggable::Static::info("Projecting reference solution on coarse mesh.");
    OGProjection<std::complex<double> > ogProjection; ogProjection.project_global(&space, &ref_sln, &sln); 
   
    // View the coarse mesh solution and polynomial orders.
    RealFilter real(&sln);
    MagFilter<double> magn(&real);
    ValFilter limited_magn(&magn, 0.0, 4e3);
    char title[100];
    sprintf(title, "Electric field, adaptivity step %d", as);
    eview.set_title(title);
    //eview.set_min_max_range(0.0, 4e3);
    eview.show(&limited_magn);
    sprintf(title, "Polynomial orders, adaptivity step %d", as);
    oview.set_title(title);
    oview.show(&space);

    // Calculate element errors and total error estimate.
    Hermes::Mixins::Loggable::Static::info("Calculating error estimate."); 
    Adapt<std::complex<double> >* adaptivity = new Adapt<std::complex<double> >(&space);

    // Set custom error form and calculate error estimate.
    CustomErrorForm cef(kappa);
    adaptivity->set_error_form(0, 0, &cef);
    double err_est_rel = adaptivity->calc_err_est(&sln, &ref_sln) * 100;

    // Report results.
    Hermes::Mixins::Loggable::Static::info("ndof_coarse: %d, ndof_fine: %d, err_est_rel: %g%%", 
      Space<std::complex<double> >::get_num_dofs(&space), 
      Space<std::complex<double> >::get_num_dofs(ref_space), err_est_rel);

    // Time measurement.
    cpu_time.tick();

    // Add entry to DOF and CPU convergence graphs.
    graph_dof.add_values(Space<std::complex<double> >::get_num_dofs(&space), err_est_rel);
    graph_dof.save("conv_dof_est.dat");
    graph_cpu.add_values(cpu_time.accumulated(), err_est_rel);
    graph_cpu.save("conv_cpu_est.dat");

    // If err_est too large, adapt the mesh.
    if (err_est_rel < ERR_STOP) done = true;
    else 
    {
      Hermes::Mixins::Loggable::Static::info("Adapting coarse mesh.");
      done = adaptivity->adapt(&selector, THRESHOLD, STRATEGY, MESH_REGULARITY);
    }
    if (space.get_num_dofs() >= NDOF_STOP) done = true;

    delete adaptivity;
    if(!done)
    {
      delete ref_space->get_mesh();
      delete ref_space;
    }
    
    // Increase counter.
    as++;
  }
  while (done == false);
  
  Hermes::Mixins::Loggable::Static::info("Total running time: %g s", cpu_time.accumulated());

  RealFilter ref_real(&sln);
  MagFilter<double> ref_magn(&ref_real);
  ValFilter ref_limited_magn(&ref_magn, 0.0, 4e3);
  eview.set_title("Fine mesh solution - magnitude");
  eview.show(&ref_limited_magn);

  // Output solution in VTK format.
  Linearizer lin;
  bool mode_3D = true;
  lin.save_solution_vtk(&ref_limited_magn, "sln.vtk", "Magnitude of E", mode_3D);
  Hermes::Mixins::Loggable::Static::info("Solution in VTK format saved to file %s.", "sln.vtk");

  // Wait for all views to be closed.
  View::wait();
  return 0;
}