<commit_before>//--------------------------------------------------------------------
//    $Id$
//
//    Copyright (C) 2005, 2006, 2010 by the deal.II authors
//
//    This file is subject to QPL and may not be  distributed
//    without copyright and license information. Please refer
//    to the file deal.II/doc/license.html for the  text  and
//    further information on this license.
//
//--------------------------------------------------------------------

// Test, whether differential operators produce a cochain complex on
// the standard Hilbert space sequence

#include "../tests.h"
#include "../lib/test_grids.h"

#include <deal.II/base/logstream.h>

#include <deal.II/lac/matrix_block.h>
#include <deal.II/lac/sparse_matrix.h>
#include <deal.II/lac/block_sparsity_pattern.h>
#include <deal.II/lac/compressed_sparsity_pattern.h>
#include <deal.II/lac/solver_cg.h>
#include <deal.II/lac/precondition.h>

#include <deal.II/dofs/dof_handler.h>
#include <deal.II/dofs/dof_renumbering.h>
#include <deal.II/dofs/dof_tools.h>

#include <deal.II/fe/fe_q.h>
#include <deal.II/fe/fe_nedelec.h>
#include <deal.II/fe/fe_raviart_thomas.h>
#include <deal.II/fe/fe_dgq.h>
#include <deal.II/fe/fe_system.h>
#include <deal.II/fe/fe_values.h>
#include <deal.II/fe/mapping_q1.h>

#include <deal.II/numerics/mesh_worker_info.h>
#include <deal.II/numerics/mesh_worker_assembler.h>
#include <deal.II/numerics/mesh_worker_loop.h>

#include <deal.II/integrators/l2.h>
#include <deal.II/integrators/divergence.h>
#include <deal.II/integrators/laplace.h>
#include <deal.II/integrators/maxwell.h>
#include <deal.II/integrators/elasticity.h>

using namespace LocalIntegrators;

const bool debugging = false;

template <int dim>
void cell_matrix(
  MeshWorker::DoFInfo<dim>& dinfo,
  typename MeshWorker::IntegrationInfo<dim,dim>& info)
{
  unsigned int dm=0; // Matrix index
  unsigned int de=0; // Element index
  
  L2::mass_matrix(dinfo.matrix(dm++,false).matrix, info.fe_values(de));
  Laplace::cell_matrix(dinfo.matrix(dm++,false).matrix, info.fe_values(de));
  Divergence::gradient_matrix(dinfo.matrix(dm++,false).matrix, info.fe_values(de), info.fe_values(de+1));

  
  ++de;
  L2::mass_matrix(dinfo.matrix(dm++,false).matrix, info.fe_values(de));
  Maxwell::curl_curl_matrix(dinfo.matrix(dm++,false).matrix, info.fe_values(de));
  Maxwell::curl_matrix(dinfo.matrix(dm++,false).matrix, info.fe_values(de), info.fe_values(de+1));

  if (dim>2)
    {
      ++de;
      L2::mass_matrix(dinfo.matrix(dm++,false).matrix, info.fe_values(de));
      Divergence::grad_div_matrix(dinfo.matrix(dm++,false).matrix, info.fe_values(de));
      Divergence::cell_matrix(dinfo.matrix(dm++,false).matrix, info.fe_values(de), info.fe_values(de+1));
    }
  
  ++de;
  L2::mass_matrix(dinfo.matrix(dm++,false).matrix, info.fe_values(de));
}



template <int dim>
void
test_cochain(const Triangulation<dim>& tr, const FiniteElement<dim>& fe)
{
  MappingQ1<dim> mapping;
				   // Initialize DofHandler for a
				   // block system with local blocks
  DoFHandler<dim> dof(tr);
  dof.distribute_dofs(fe);
  DoFRenumbering::component_wise(dof);
  dof.initialize_local_block_info();

				   // Initialize hanging node constraints
  ConstraintMatrix constraints;
//  DoFTools::make_zero_boundary_constraints(dof, constraints);
  DoFTools::make_hanging_node_constraints (dof, constraints);
  constraints.close();
				   // Setup sparsity pattern for the
				   // whole system
  BlockSparsityPattern sparsity;
  sparsity.reinit(dof.block_info().global().size(),
		  dof.block_info().global().size());
  BlockCompressedSparsityPattern c_sparsity(dof.block_info().global(),
					    dof.block_info().global());
  DoFTools::make_flux_sparsity_pattern(dof, c_sparsity, constraints, false);
  sparsity.copy_from(c_sparsity);

				   // Setup matrices
  MatrixBlockVector<SparseMatrix<double> > matrices;

  unsigned int d=0;
  matrices.add(d,d,"mass-H1");
  matrices.add(d,d,"Laplacian");
  matrices.add(d+1,d,"grad");
  
  ++d;
  matrices.add(d,d,"mass-Hcurl");
  matrices.add(d,d,"curl-curl");
  matrices.add(d+1,d,"curl");
  
  if (dim>2)
    {
      ++d;
      matrices.add(d,d,"mass-Hdiv");
      matrices.add(d,d,"grad-div");
      matrices.add(d+1,d,"div");
    }

  ++d;
  matrices.add(d,d,"mass-L2");

  matrices.reinit(sparsity);

  if (debugging)
    for (unsigned int i=0;i<matrices.size();++i)
      deallog << "Block " << '(' << matrices.block(i).row
	      << ',' << matrices.block(i).column << ") "
	      << std::setw(3) << std::right << matrices.matrix(i).m() << std::left
	      << 'x' << std::setw(3) << matrices.matrix(i).n()
	      << ' ' << matrices.name(i) << std::endl;
  
				   // Build matrices
  
  MeshWorker::IntegrationInfoBox<dim> info_box;
  UpdateFlags update_flags = update_values | update_gradients;
  info_box.add_update_flags_cell(update_flags);
  info_box.initialize(fe, mapping, &dof.block_info());
  if (debugging)
    deallog << "Infobox ready" << std::endl;
  
  MeshWorker::DoFInfo<dim> dof_info(dof.block_info());
  
  MeshWorker::Assembler::MatrixLocalBlocksToGlobalBlocks<SparseMatrix<double> > assembler;
  assembler.initialize(&dof.block_info(), matrices);
  assembler.initialize(constraints);
  
  MeshWorker::integration_loop<dim, dim>(
    dof.begin_active(), dof.end(), dof_info, info_box,
    cell_matrix<dim>, 0, 0, assembler);
  
  for (unsigned int b=0;b<matrices.size();++b)
    if (b%3 == 0)
      for (unsigned int i=0;i<matrices.block(b).matrix.m();++i)
	if (matrices.block(b).matrix.diag_element(i) == 0.)
	  matrices.block(b).matrix.diag_element(i) = 1.;
  if (debugging)
    deallog << "Matrices ready" << std::endl;
				   // Set up vectors
  BlockVector<double> source(dof.block_info().global());
  BlockVector<double> result1(dof.block_info().global());
  BlockVector<double> result2(dof.block_info().global());
  BlockVector<double> aux(dof.block_info().global());
  for (unsigned int i=0;i<source.size();++i)
    source(i) = i%5;
<commit_msg>avoid deprecations<commit_after>//--------------------------------------------------------------------
//    $Id$
//
//    Copyright (C) 2005, 2006, 2010, 2013 by the deal.II authors
//
//    This file is subject to QPL and may not be  distributed
//    without copyright and license information. Please refer
//    to the file deal.II/doc/license.html for the  text  and
//    further information on this license.
//
//--------------------------------------------------------------------

// Test, whether differential operators produce a cochain complex on
// the standard Hilbert space sequence

#include "../tests.h"
#include "../lib/test_grids.h"

#include <deal.II/base/logstream.h>

#include <deal.II/lac/matrix_block.h>
#include <deal.II/lac/sparse_matrix.h>
#include <deal.II/lac/block_sparsity_pattern.h>
#include <deal.II/lac/compressed_sparsity_pattern.h>
#include <deal.II/lac/solver_cg.h>
#include <deal.II/lac/precondition.h>

#include <deal.II/dofs/dof_handler.h>
#include <deal.II/dofs/dof_renumbering.h>
#include <deal.II/dofs/dof_tools.h>

#include <deal.II/fe/fe_q.h>
#include <deal.II/fe/fe_nedelec.h>
#include <deal.II/fe/fe_raviart_thomas.h>
#include <deal.II/fe/fe_dgq.h>
#include <deal.II/fe/fe_system.h>
#include <deal.II/fe/fe_values.h>
#include <deal.II/fe/mapping_q1.h>

#include <deal.II/meshworker/integration_info.h>
#include <deal.II/meshworker/assembler.h>
#include <deal.II/meshworker/loop.h>

#include <deal.II/integrators/l2.h>
#include <deal.II/integrators/divergence.h>
#include <deal.II/integrators/laplace.h>
#include <deal.II/integrators/maxwell.h>
#include <deal.II/integrators/elasticity.h>

using namespace LocalIntegrators;

const bool debugging = false;

template <int dim>
void cell_matrix(
  MeshWorker::DoFInfo<dim>& dinfo,
  typename MeshWorker::IntegrationInfo<dim,dim>& info)
{
  unsigned int dm=0; // Matrix index
  unsigned int de=0; // Element index
  
  L2::mass_matrix(dinfo.matrix(dm++,false).matrix, info.fe_values(de));
  Laplace::cell_matrix(dinfo.matrix(dm++,false).matrix, info.fe_values(de));
  Divergence::gradient_matrix(dinfo.matrix(dm++,false).matrix, info.fe_values(de), info.fe_values(de+1));

  
  ++de;
  L2::mass_matrix(dinfo.matrix(dm++,false).matrix, info.fe_values(de));
  Maxwell::curl_curl_matrix(dinfo.matrix(dm++,false).matrix, info.fe_values(de));
  Maxwell::curl_matrix(dinfo.matrix(dm++,false).matrix, info.fe_values(de), info.fe_values(de+1));

  if (dim>2)
    {
      ++de;
      L2::mass_matrix(dinfo.matrix(dm++,false).matrix, info.fe_values(de));
      Divergence::grad_div_matrix(dinfo.matrix(dm++,false).matrix, info.fe_values(de));
      Divergence::cell_matrix(dinfo.matrix(dm++,false).matrix, info.fe_values(de), info.fe_values(de+1));
    }
  
  ++de;
  L2::mass_matrix(dinfo.matrix(dm++,false).matrix, info.fe_values(de));
}



template <int dim>
void
test_cochain(const Triangulation<dim>& tr, const FiniteElement<dim>& fe)
{
  MappingQ1<dim> mapping;
				   // Initialize DofHandler for a
				   // block system with local blocks
  DoFHandler<dim> dof(tr);
  dof.distribute_dofs(fe);
  DoFRenumbering::component_wise(dof);
  dof.initialize_local_block_info();

				   // Initialize hanging node constraints
  ConstraintMatrix constraints;
//  DoFTools::make_zero_boundary_constraints(dof, constraints);
  DoFTools::make_hanging_node_constraints (dof, constraints);
  constraints.close();
				   // Setup sparsity pattern for the
				   // whole system
  BlockSparsityPattern sparsity;
  sparsity.reinit(dof.block_info().global().size(),
		  dof.block_info().global().size());
  BlockCompressedSparsityPattern c_sparsity(dof.block_info().global(),
					    dof.block_info().global());
  DoFTools::make_flux_sparsity_pattern(dof, c_sparsity, constraints, false);
  sparsity.copy_from(c_sparsity);

				   // Setup matrices
  MatrixBlockVector<SparseMatrix<double> > matrices;

  unsigned int d=0;
  matrices.add(d,d,"mass-H1");
  matrices.add(d,d,"Laplacian");
  matrices.add(d+1,d,"grad");
  
  ++d;
  matrices.add(d,d,"mass-Hcurl");
  matrices.add(d,d,"curl-curl");
  matrices.add(d+1,d,"curl");
  
  if (dim>2)
    {
      ++d;
      matrices.add(d,d,"mass-Hdiv");
      matrices.add(d,d,"grad-div");
      matrices.add(d+1,d,"div");
    }

  ++d;
  matrices.add(d,d,"mass-L2");

  matrices.reinit(sparsity);

  if (debugging)
    for (unsigned int i=0;i<matrices.size();++i)
      deallog << "Block " << '(' << matrices.block(i).row
	      << ',' << matrices.block(i).column << ") "
	      << std::setw(3) << std::right << matrices.matrix(i).m() << std::left
	      << 'x' << std::setw(3) << matrices.matrix(i).n()
	      << ' ' << matrices.name(i) << std::endl;
  
				   // Build matrices
  
  MeshWorker::IntegrationInfoBox<dim> info_box;
  UpdateFlags update_flags = update_values | update_gradients;
  info_box.add_update_flags_cell(update_flags);
  info_box.initialize(fe, mapping, &dof.block_info());
  if (debugging)
    deallog << "Infobox ready" << std::endl;
  
  MeshWorker::DoFInfo<dim> dof_info(dof.block_info());
  
  MeshWorker::Assembler::MatrixLocalBlocksToGlobalBlocks<SparseMatrix<double> > assembler;
  assembler.initialize(&dof.block_info(), matrices);
  assembler.initialize(constraints);
  
  MeshWorker::loop<dim, dim, MeshWorker::DoFInfo<dim>, MeshWorker::IntegrationInfoBox<dim> >(
    dof.begin_active(), dof.end(), dof_info, info_box,
    cell_matrix<dim>, 0, 0, assembler);
  
  for (unsigned int b=0;b<matrices.size();++b)
    if (b%3 == 0)
      for (unsigned int i=0;i<matrices.block(b).matrix.m();++i)
	if (matrices.block(b).matrix.diag_element(i) == 0.)
	  matrices.block(b).matrix.diag_element(i) = 1.;
  if (debugging)
    deallog << "Matrices ready" << std::endl;
				   // Set up vectors
  BlockVector<double> source(dof.block_info().global());
  BlockVector<double> result1(dof.block_info().global());
  BlockVector<double> result2(dof.block_info().global());
  BlockVector<double> aux(dof.block_info().global());
  for (unsigned int i=0;i<source.size();++i)
    source(i) = i%5;

{'fraction_non_alphanumeric': 0.1894640326748354, 'fraction_numerical': 0.008168708843877636, 'mean_word_length': 3.9578512396694214, 'pattern_counts': {'":': 0, '<': 183, '<?xml version=': 0, '>': 99, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 0, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'dir': 'git-commits-cleaned', 'id': '6321954', 'n_tokens_mistral': 4616, 'n_tokens_neox': 4349, 'n_words': 680}