Hi I'm studying step-20 tutorial example, watching Wolfgang's Lectures
In step20, I would like to solve that Eq in the 3D mesh that I made by
myself(I attached B_coord3D.inp that is UCD format)
But, It didn't work out.
The error message is like this....
[ 50%] Built target step-20
make[3]: Warning: File `step-20' has modification time 81 s in the future
[100%] Run step-20 with Debug configuration
Number of active cells: 36864
Total number of cells: 42048
Number of degrees of freedom: 156672 (119808+36864)
53 CG Schur complement iterations to obtain convergence.
--------------------------------------------------------
An error occurred in line <790> of file
</user2/hanks318/dealii/dealii-8.3.0/source/fe/mapping_q1.cc> in function
void dealii::MappingQ1<dim, spacedim>::fill_fe_values(const typename
dealii::Triangulation<dim, spacedim>::cell_iterator&, const
dealii::Quadrature<dim>&, typename dealii::Mapping<dim,
spacedim>::InternalDataBase&, std::vector<dealii::Point<spacedim, double>,
std::allocator<dealii::Point<spacedim, double> > >&, std::vector<double,
std::allocator<double> >&, std::vector<dealii::DerivativeForm<1, spacedim,
dim, double>, std::allocator<dealii::DerivativeForm<1, spacedim, dim,
double> > >&, std::vector<dealii::DerivativeForm<2, dim, spacedim, double>,
std::allocator<dealii::DerivativeForm<2, dim, spacedim, double> > >&,
std::vector<dealii::DerivativeForm<1, spacedim, dim, double>,
std::allocator<dealii::DerivativeForm<1, spacedim, dim, double> > >&,
std::vector<dealii::Point<spacedim, double>,
std::allocator<dealii::Point<spacedim, double> > >&,
dealii::CellSimilarity::Similarity&) const [with int dim = 3, int spacedim
= 3]
The violated condition was:
det > 1e-12*Utilities::fixed_power<dim>(cell->diameter()/
std::sqrt(double(dim)))
The name and call sequence of the exception was:
(typename Mapping<dim,spacedim>::ExcDistortedMappedCell(cell->center(),
det, point))
Additional Information:
The image of the mapping applied to cell with center [-1.61075 -0.761538
0.649647] is distorted. The cell geometry or the mapping are invalid,
giving a non-positive volume fraction of -1.49556e-07 in quadrature point 5.
Stacktrace:
-----------
#0 /user2/hanks318/dealii/dealii_pre/lib/libdeal_II.g.so.8.3.0:
dealii::MappingQ1<3,
3>::fill_fe_values(dealii::TriaIterator<dealii::CellAccessor<3, 3> >
const&, dealii::Quadrature<3> const&, dealii::Mapping<3,
3>::InternalDataBase&, std::vector<dealii::Point<3, double>,
std::allocator<dealii::Point<3, double> > >&, std::vector<double,
std::allocator<double> >&, std::vector<dealii::DerivativeForm<1, 3, 3,
double>, std::allocator<dealii::DerivativeForm<1, 3, 3, double> > >&,
std::vector<dealii::DerivativeForm<2, 3, 3, double>,
std::allocator<dealii::DerivativeForm<2, 3, 3, double> > >&,
std::vector<dealii::DerivativeForm<1, 3, 3, double>,
std::allocator<dealii::DerivativeForm<1, 3, 3, double> > >&,
std::vector<dealii::Point<3, double>, std::allocator<dealii::Point<3,
double> > >&, dealii::CellSimilarity::Similarity&) const
#1 /user2/hanks318/dealii/dealii_pre/lib/libdeal_II.g.so.8.3.0:
dealii::FEValues<3, 3>::do_reinit()
#2 /user2/hanks318/dealii/dealii_pre/lib/libdeal_II.g.so.8.3.0: void
dealii::FEValues<3, 3>::reinit<dealii::DoFHandler<3, 3>,
false>(dealii::TriaIterator<dealii::DoFCellAccessor<dealii::DoFHandler<3,
3>, false> > const&)
#3 /user2/hanks318/dealii/dealii_pre/lib/libdeal_II.g.so.8.3.0: void
dealii::hp::FEValues<3, 3>::reinit<dealii::DoFHandler<3, 3>,
false>(dealii::TriaIterator<dealii::DoFCellAccessor<dealii::DoFHandler<3,
3>, false> >, unsigned int, unsigned int, unsigned int)
#4 /user2/hanks318/dealii/dealii_pre/lib/libdeal_II.g.so.8.3.0: void
dealii::internal::DataOut::ParallelDataBase<3,
3>::reinit_all_fe_values<dealii::DoFHandler<3, 3>
>(std::vector<boost::shared_ptr<dealii::internal::DataOut::DataEntryBase<dealii::DoFHandler<3,
>
3> > >,
std::allocator<boost::shared_ptr<dealii::internal::DataOut::DataEntryBase<dealii::DoFHandler<3,
3> > > > >&, dealii::TriaIterator<dealii::CellAccessor<3, 3> > const&,
unsigned int)
#5 /user2/hanks318/dealii/dealii_pre/lib/libdeal_II.g.so.8.3.0:
dealii::DataOut<3, dealii::DoFHandler<3, 3>
>::build_one_patch(std::pair<dealii::TriaIterator<dealii::CellAccessor<3,
3> >, unsigned int> const*, dealii::internal::DataOut::ParallelData<3, 3>&,
dealii::DataOutBase::Patch<3, 3>&, dealii::DataOut<3, dealii::DoFHandler<3,
3> >::CurvedCellRegion, std::vector<dealii::DataOutBase::Patch<3, 3>,
std::allocator<dealii::DataOutBase::Patch<3, 3> > >&)
#6 /user2/hanks318/dealii/dealii_pre/lib/libdeal_II.g.so.8.3.0:
boost::_mfi::mf5<void, dealii::DataOut<3, dealii::DoFHandler<3, 3> >,
std::pair<dealii::TriaIterator<dealii::CellAccessor<3, 3> >, unsigned int>
const*, dealii::internal::DataOut::ParallelData<3, 3>&,
dealii::DataOutBase::Patch<3, 3>&, dealii::DataOut<3, dealii::DoFHandler<3,
3> >::CurvedCellRegion, std::vector<dealii::DataOutBase::Patch<3, 3>,
std::allocator<dealii::DataOutBase::Patch<3, 3> >
>&>::operator()(dealii::DataOut<3, dealii::DoFHandler<3, 3> >*,
std::pair<dealii::TriaIterator<dealii::CellAccessor<3, 3> >, unsigned int>
const*, dealii::internal::DataOut::ParallelData<3, 3>&,
dealii::DataOutBase::Patch<3, 3>&, dealii::DataOut<3, dealii::DoFHandler<3,
3> >::CurvedCellRegion, std::vector<dealii::DataOutBase::Patch<3, 3>,
std::allocator<dealii::DataOutBase::Patch<3, 3> > >&) const
#7 /user2/hanks318/dealii/dealii_pre/lib/libdeal_II.g.so.8.3.0: void
boost::_bi::list6<boost::_bi::value<dealii::DataOut<3,
dealii::DoFHandler<3, 3> >*>, boost::arg<1>, boost::arg<2>, boost::arg<3>,
boost::_bi::value<dealii::DataOut<3, dealii::DoFHandler<3, 3>
>::CurvedCellRegion>,
boost::reference_wrapper<std::vector<dealii::DataOutBase::Patch<3, 3>,
std::allocator<dealii::DataOutBase::Patch<3, 3> > > >
>::operator()<boost::_mfi::mf5<void, dealii::DataOut<3,
dealii::DoFHandler<3, 3> >,
std::pair<dealii::TriaIterator<dealii::CellAccessor<3, 3> >, unsigned int>
const*, dealii::internal::DataOut::ParallelData<3, 3>&,
dealii::DataOutBase::Patch<3, 3>&, dealii::DataOut<3, dealii::DoFHandler<3,
3> >::CurvedCellRegion, std::vector<dealii::DataOutBase::Patch<3, 3>,
std::allocator<dealii::DataOutBase::Patch<3, 3> > >&>,
boost::_bi::list3<std::pair<dealii::TriaIterator<dealii::CellAccessor<3, 3>
>, unsigned int>* const&, dealii::internal::DataOut::ParallelData<3, 3>&,
dealii::DataOutBase::Patch<3, 3>&> >(boost::_bi::type<void>,
boost::_mfi::mf5<void, dealii::DataOut<3, dealii::DoFHandler<3, 3> >,
std::pair<dealii::TriaIterator<dealii::CellAccessor<3, 3> >, unsigned int>
const*, dealii::internal::DataOut::ParallelData<3, 3>&,
dealii::DataOutBase::Patch<3, 3>&, dealii::DataOut<3, dealii::DoFHandler<3,
3> >::CurvedCellRegion, std::vector<dealii::DataOutBase::Patch<3, 3>,
std::allocator<dealii::DataOutBase::Patch<3, 3> > >&>&,
boost::_bi::list3<std::pair<dealii::TriaIterator<dealii::CellAccessor<3, 3>
>, unsigned int>* const&, dealii::internal::DataOut::ParallelData<3, 3>&,
dealii::DataOutBase::Patch<3, 3>&>&, int)
#8 /user2/hanks318/dealii/dealii_pre/lib/libdeal_II.g.so.8.3.0: void
boost::_bi::bind_t<void, boost::_mfi::mf5<void, dealii::DataOut<3,
dealii::DoFHandler<3, 3> >,
std::pair<dealii::TriaIterator<dealii::CellAccessor<3, 3> >, unsigned int>
const*, dealii::internal::DataOut::ParallelData<3, 3>&,
dealii::DataOutBase::Patch<3, 3>&, dealii::DataOut<3, dealii::DoFHandler<3,
3> >::CurvedCellRegion, std::vector<dealii::DataOutBase::Patch<3, 3>,
std::allocator<dealii::DataOutBase::Patch<3, 3> > >&>,
boost::_bi::list6<boost::_bi::value<dealii::DataOut<3,
dealii::DoFHandler<3, 3> >*>, boost::arg<1>, boost::arg<2>, boost::arg<3>,
boost::_bi::value<dealii::DataOut<3, dealii::DoFHandler<3, 3>
>::CurvedCellRegion>,
boost::reference_wrapper<std::vector<dealii::DataOutBase::Patch<3, 3>,
std::allocator<dealii::DataOutBase::Patch<3, 3> > > > >
>::operator()<std::pair<dealii::TriaIterator<dealii::CellAccessor<3, 3> >,
unsigned int>* const, dealii::internal::DataOut::ParallelData<3, 3>,
dealii::DataOutBase::Patch<3, 3>
>(std::pair<dealii::TriaIterator<dealii::CellAccessor<3, 3> >, unsigned
int>* const&, dealii::internal::DataOut::ParallelData<3, 3>&,
dealii::DataOutBase::Patch<3, 3>&)
#9 /user2/hanks318/dealii/dealii_pre/lib/libdeal_II.g.so.8.3.0:
boost::detail::function::void_function_obj_invoker3<boost::_bi::bind_t<void,
boost::_mfi::mf5<void, dealii::DataOut<3, dealii::DoFHandler<3, 3> >,
std::pair<dealii::TriaIterator<dealii::CellAccessor<3, 3> >, unsigned int>
const*, dealii::internal::DataOut::ParallelData<3, 3>&,
dealii::DataOutBase::Patch<3, 3>&, dealii::DataOut<3, dealii::DoFHandler<3,
3> >::CurvedCellRegion, std::vector<dealii::DataOutBase::Patch<3, 3>,
std::allocator<dealii::DataOutBase::Patch<3, 3> > >&>,
boost::_bi::list6<boost::_bi::value<dealii::DataOut<3,
dealii::DoFHandler<3, 3> >*>, boost::arg<1>, boost::arg<2>, boost::arg<3>,
boost::_bi::value<dealii::DataOut<3, dealii::DoFHandler<3, 3>
>::CurvedCellRegion>,
boost::reference_wrapper<std::vector<dealii::DataOutBase::Patch<3, 3>,
std::allocator<dealii::DataOutBase::Patch<3, 3> > > > > >, void,
std::pair<dealii::TriaIterator<dealii::CellAccessor<3, 3> >, unsigned int>*
const&, dealii::internal::DataOut::ParallelData<3, 3>&,
dealii::DataOutBase::Patch<3,
3>&>::invoke(boost::detail::function::function_buffer&,
std::pair<dealii::TriaIterator<dealii::CellAccessor<3, 3> >, unsigned int>*
const&, dealii::internal::DataOut::ParallelData<3, 3>&,
dealii::DataOutBase::Patch<3, 3>&)
#10 /user2/hanks318/dealii/dealii_pre/lib/libdeal_II.g.so.8.3.0:
boost::function3<void,
std::pair<dealii::TriaIterator<dealii::CellAccessor<3, 3> >, unsigned int>*
const&, dealii::internal::DataOut::ParallelData<3, 3>&,
dealii::DataOutBase::Patch<3,
3>&>::operator()(std::pair<dealii::TriaIterator<dealii::CellAccessor<3, 3>
>, unsigned int>* const&, dealii::internal::DataOut::ParallelData<3, 3>&,
dealii::DataOutBase::Patch<3, 3>&) const
#11 /user2/hanks318/dealii/dealii_pre/lib/libdeal_II.g.so.8.3.0:
dealii::WorkStream::internal::Implementation3::WorkerAndCopier<std::pair<dealii::TriaIterator<dealii::CellAccessor<3,
3> >, unsigned int>*, dealii::internal::DataOut::ParallelData<3, 3>,
dealii::DataOutBase::Patch<3, 3>
>::operator()(tbb::blocked_range<__gnu_cxx::__normal_iterator<std::pair<dealii::TriaIterator<dealii::CellAccessor<3,
>
3> >, unsigned int>* const*,
std::vector<std::pair<dealii::TriaIterator<dealii::CellAccessor<3, 3> >,
unsigned int>*,
std::allocator<std::pair<dealii::TriaIterator<dealii::CellAccessor<3, 3> >,
unsigned int>*> > > > const&)
#12 /user2/hanks318/dealii/dealii_pre/lib/libdeal_II.g.so.8.3.0:
boost::_mfi::mf1<void,
dealii::WorkStream::internal::Implementation3::WorkerAndCopier<std::pair<dealii::TriaIterator<dealii::CellAccessor<3,
3> >, unsigned int>*, dealii::internal::DataOut::ParallelData<3, 3>,
dealii::DataOutBase::Patch<3, 3> >,
tbb::blocked_range<__gnu_cxx::__normal_iterator<std::pair<dealii::TriaIterator<dealii::CellAccessor<3,
3> >, unsigned int>* const*,
std::vector<std::pair<dealii::TriaIterator<dealii::CellAccessor<3, 3> >,
unsigned int>*,
std::allocator<std::pair<dealii::TriaIterator<dealii::CellAccessor<3, 3> >,
unsigned int>*> > > >
const&>::operator()(dealii::WorkStream::internal::Implementation3::WorkerAndCopier<std::pair<dealii::TriaIterator<dealii::CellAccessor<3,
3> >, unsigned int>*, dealii::internal::DataOut::ParallelData<3, 3>,
dealii::DataOutBase::Patch<3, 3> >&,
tbb::blocked_range<__gnu_cxx::__normal_iterator<std::pair<dealii::TriaIterator<dealii::CellAccessor<3,
3> >, unsigned int>* const*,
std::vector<std::pair<dealii::TriaIterator<dealii::CellAccessor<3, 3> >,
unsigned int>*,
std::allocator<std::pair<dealii::TriaIterator<dealii::CellAccessor<3, 3> >,
unsigned int>*> > > > const&) const
#13 /user2/hanks318/dealii/dealii_pre/lib/libdeal_II.g.so.8.3.0: void
boost::_bi::list2<boost::reference_wrapper<dealii::WorkStream::internal::Implementation3::WorkerAndCopier<std::pair<dealii::TriaIterator<dealii::CellAccessor<3,
3> >, unsigned int>*, dealii::internal::DataOut::ParallelData<3, 3>,
dealii::DataOutBase::Patch<3, 3> > >, boost::arg<1>
>::operator()<boost::_mfi::mf1<void,
dealii::WorkStream::internal::Implementation3::WorkerAndCopier<std::pair<dealii::TriaIterator<dealii::CellAccessor<3,
3> >, unsigned int>*, dealii::internal::DataOut::ParallelData<3, 3>,
dealii::DataOutBase::Patch<3, 3> >,
tbb::blocked_range<__gnu_cxx::__normal_iterator<std::pair<dealii::TriaIterator<dealii::CellAccessor<3,
3> >, unsigned int>* const*,
std::vector<std::pair<dealii::TriaIterator<dealii::CellAccessor<3, 3> >,
unsigned int>*,
std::allocator<std::pair<dealii::TriaIterator<dealii::CellAccessor<3, 3> >,
unsigned int>*> > > > const&>,
boost::_bi::list1<tbb::blocked_range<__gnu_cxx::__normal_iterator<std::pair<dealii::TriaIterator<dealii::CellAccessor<3,
3> >, unsigned int>* const*,
std::vector<std::pair<dealii::TriaIterator<dealii::CellAccessor<3, 3> >,
unsigned int>*,
std::allocator<std::pair<dealii::TriaIterator<dealii::CellAccessor<3, 3> >,
unsigned int>*> > > >&> >(boost::_bi::type<void>, boost::_mfi::mf1<void,
dealii::WorkStream::internal::Implementation3::WorkerAndCopier<std::pair<dealii::TriaIterator<dealii::CellAccessor<3,
3> >, unsigned int>*, dealii::internal::DataOut::ParallelData<3, 3>,
dealii::DataOutBase::Patch<3, 3> >,
tbb::blocked_range<__gnu_cxx::__normal_iterator<std::pair<dealii::TriaIterator<dealii::CellAccessor<3,
3> >, unsigned int>* const*,
std::vector<std::pair<dealii::TriaIterator<dealii::CellAccessor<3, 3> >,
unsigned int>*,
std::allocator<std::pair<dealii::TriaIterator<dealii::CellAccessor<3, 3> >,
unsigned int>*> > > > const&> const&,
boost::_bi::list1<tbb::blocked_range<__gnu_cxx::__normal_iterator<std::pair<dealii::TriaIterator<dealii::CellAccessor<3,
3> >, unsigned int>* const*,
std::vector<std::pair<dealii::TriaIterator<dealii::CellAccessor<3, 3> >,
unsigned int>*,
std::allocator<std::pair<dealii::TriaIterator<dealii::CellAccessor<3, 3> >,
unsigned int>*> > > >&>&, int) const
#14 /user2/hanks318/dealii/dealii_pre/lib/libdeal_II.g.so.8.3.0: void
boost::_bi::bind_t<void, boost::_mfi::mf1<void,
dealii::WorkStream::internal::Implementation3::WorkerAndCopier<std::pair<dealii::TriaIterator<dealii::CellAccessor<3,
3> >, unsigned int>*, dealii::internal::DataOut::ParallelData<3, 3>,
dealii::DataOutBase::Patch<3, 3> >,
tbb::blocked_range<__gnu_cxx::__normal_iterator<std::pair<dealii::TriaIterator<dealii::CellAccessor<3,
3> >, unsigned int>* const*,
std::vector<std::pair<dealii::TriaIterator<dealii::CellAccessor<3, 3> >,
unsigned int>*,
std::allocator<std::pair<dealii::TriaIterator<dealii::CellAccessor<3, 3> >,
unsigned int>*> > > > const&>,
boost::_bi::list2<boost::reference_wrapper<dealii::WorkStream::internal::Implementation3::WorkerAndCopier<std::pair<dealii::TriaIterator<dealii::CellAccessor<3,
3> >, unsigned int>*, dealii::internal::DataOut::ParallelData<3, 3>,
dealii::DataOutBase::Patch<3, 3> > >, boost::arg<1> >
>::operator()<tbb::blocked_range<__gnu_cxx::__normal_iterator<std::pair<dealii::TriaIterator<dealii::CellAccessor<3,
>
3> >, unsigned int>* const*,
std::vector<std::pair<dealii::TriaIterator<dealii::CellAccessor<3, 3> >,
unsigned int>*,
std::allocator<std::pair<dealii::TriaIterator<dealii::CellAccessor<3, 3> >,
unsigned int>*> > > >
>(tbb::blocked_range<__gnu_cxx::__normal_iterator<std::pair<dealii::TriaIterator<dealii::CellAccessor<3,
>
3> >, unsigned int>* const*,
std::vector<std::pair<dealii::TriaIterator<dealii::CellAccessor<3, 3> >,
unsigned int>*,
std::allocator<std::pair<dealii::TriaIterator<dealii::CellAccessor<3, 3> >,
unsigned int>*> > > >&) const
#15 /user2/hanks318/dealii/dealii_pre/lib/libdeal_II.g.so.8.3.0:
tbb::interface6::internal::start_for<tbb::blocked_range<__gnu_cxx::__normal_iterator<std::pair<dealii::TriaIterator<dealii::CellAccessor<3,
3> >, unsigned int>* const*,
std::vector<std::pair<dealii::TriaIterator<dealii::CellAccessor<3, 3> >,
unsigned int>*,
std::allocator<std::pair<dealii::TriaIterator<dealii::CellAccessor<3, 3> >,
unsigned int>*> > > >, boost::_bi::bind_t<void, boost::_mfi::mf1<void,
dealii::WorkStream::internal::Implementation3::WorkerAndCopier<std::pair<dealii::TriaIterator<dealii::CellAccessor<3,
3> >, unsigned int>*, dealii::internal::DataOut::ParallelData<3, 3>,
dealii::DataOutBase::Patch<3, 3> >,
tbb::blocked_range<__gnu_cxx::__normal_iterator<std::pair<dealii::TriaIterator<dealii::CellAccessor<3,
3> >, unsigned int>* const*,
std::vector<std::pair<dealii::TriaIterator<dealii::CellAccessor<3, 3> >,
unsigned int>*,
std::allocator<std::pair<dealii::TriaIterator<dealii::CellAccessor<3, 3> >,
unsigned int>*> > > > const&>,
boost::_bi::list2<boost::reference_wrapper<dealii::WorkStream::internal::Implementation3::WorkerAndCopier<std::pair<dealii::TriaIterator<dealii::CellAccessor<3,
3> >, unsigned int>*, dealii::internal::DataOut::ParallelData<3, 3>,
dealii::DataOutBase::Patch<3, 3> > >, boost::arg<1> > >,
tbb::auto_partitioner
const>::run_body(tbb::blocked_range<__gnu_cxx::__normal_iterator<std::pair<dealii::TriaIterator<dealii::CellAccessor<3,
3> >, unsigned int>* const*,
std::vector<std::pair<dealii::TriaIterator<dealii::CellAccessor<3, 3> >,
unsigned int>*,
std::allocator<std::pair<dealii::TriaIterator<dealii::CellAccessor<3, 3> >,
unsigned int>*> > > >&)
#16 /user2/hanks318/dealii/dealii_pre/lib/libdeal_II.g.so.8.3.0: void
tbb::interface6::internal::partition_type_base<tbb::interface6::internal::auto_partition_type>::execute<tbb::interface6::internal::start_for<tbb::blocked_range<__gnu_cxx::__normal_iterator<std::pair<dealii::TriaIterator<dealii::CellAccessor<3,
3> >, unsigned int>* const*,
std::vector<std::pair<dealii::TriaIterator<dealii::CellAccessor<3, 3> >,
unsigned int>*,
std::allocator<std::pair<dealii::TriaIterator<dealii::CellAccessor<3, 3> >,
unsigned int>*> > > >, boost::_bi::bind_t<void, boost::_mfi::mf1<void,
dealii::WorkStream::internal::Implementation3::WorkerAndCopier<std::pair<dealii::TriaIterator<dealii::CellAccessor<3,
3> >, unsigned int>*, dealii::internal::DataOut::ParallelData<3, 3>,
dealii::DataOutBase::Patch<3, 3> >,
tbb::blocked_range<__gnu_cxx::__normal_iterator<std::pair<dealii::TriaIterator<dealii::CellAccessor<3,
3> >, unsigned int>* const*,
std::vector<std::pair<dealii::TriaIterator<dealii::CellAccessor<3, 3> >,
unsigned int>*,
std::allocator<std::pair<dealii::TriaIterator<dealii::CellAccessor<3, 3> >,
unsigned int>*> > > > const&>,
boost::_bi::list2<boost::reference_wrapper<dealii::WorkStream::internal::Implementation3::WorkerAndCopier<std::pair<dealii::TriaIterator<dealii::CellAccessor<3,
3> >, unsigned int>*, dealii::internal::DataOut::ParallelData<3, 3>,
dealii::DataOutBase::Patch<3, 3> > >, boost::arg<1> > >,
tbb::auto_partitioner const>,
tbb::blocked_range<__gnu_cxx::__normal_iterator<std::pair<dealii::TriaIterator<dealii::CellAccessor<3,
3> >, unsigned int>* const*,
std::vector<std::pair<dealii::TriaIterator<dealii::CellAccessor<3, 3> >,
unsigned int>*,
std::allocator<std::pair<dealii::TriaIterator<dealii::CellAccessor<3, 3> >,
unsigned int>*> > > >
>(tbb::interface6::internal::start_for<tbb::blocked_range<__gnu_cxx::__normal_iterator<std::pair<dealii::TriaIterator<dealii::CellAccessor<3,
>
3> >, unsigned int>* const*,
std::vector<std::pair<dealii::TriaIterator<dealii::CellAccessor<3, 3> >,
unsigned int>*,
std::allocator<std::pair<dealii::TriaIterator<dealii::CellAccessor<3, 3> >,
unsigned int>*> > > >, boost::_bi::bind_t<void, boost::_mfi::mf1<void,
dealii::WorkStream::internal::Implementation3::WorkerAndCopier<std::pair<dealii::TriaIterator<dealii::CellAccessor<3,
3> >, unsigned int>*, dealii::internal::DataOut::ParallelData<3, 3>,
dealii::DataOutBase::Patch<3, 3> >,
tbb::blocked_range<__gnu_cxx::__normal_iterator<std::pair<dealii::TriaIterator<dealii::CellAccessor<3,
3> >, unsigned int>* const*,
std::vector<std::pair<dealii::TriaIterator<dealii::CellAccessor<3, 3> >,
unsigned int>*,
std::allocator<std::pair<dealii::TriaIterator<dealii::CellAccessor<3, 3> >,
unsigned int>*> > > > const&>,
boost::_bi::list2<boost::reference_wrapper<dealii::WorkStream::internal::Implementation3::WorkerAndCopier<std::pair<dealii::TriaIterator<dealii::CellAccessor<3,
3> >, unsigned int>*, dealii::internal::DataOut::ParallelData<3, 3>,
dealii::DataOutBase::Patch<3, 3> > >, boost::arg<1> > >,
tbb::auto_partitioner const>&,
tbb::blocked_range<__gnu_cxx::__normal_iterator<std::pair<dealii::TriaIterator<dealii::CellAccessor<3,
3> >, unsigned int>* const*,
std::vector<std::pair<dealii::TriaIterator<dealii::CellAccessor<3, 3> >,
unsigned int>*,
std::allocator<std::pair<dealii::TriaIterator<dealii::CellAccessor<3, 3> >,
unsigned int>*> > > >&)
#17 /user2/hanks318/dealii/dealii_pre/lib/libdeal_II.g.so.8.3.0:
tbb::interface6::internal::start_for<tbb::blocked_range<__gnu_cxx::__normal_iterator<std::pair<dealii::TriaIterator<dealii::CellAccessor<3,
3> >, unsigned int>* const*,
std::vector<std::pair<dealii::TriaIterator<dealii::CellAccessor<3, 3> >,
unsigned int>*,
std::allocator<std::pair<dealii::TriaIterator<dealii::CellAccessor<3, 3> >,
unsigned int>*> > > >, boost::_bi::bind_t<void, boost::_mfi::mf1<void,
dealii::WorkStream::internal::Implementation3::WorkerAndCopier<std::pair<dealii::TriaIterator<dealii::CellAccessor<3,
3> >, unsigned int>*, dealii::internal::DataOut::ParallelData<3, 3>,
dealii::DataOutBase::Patch<3, 3> >,
tbb::blocked_range<__gnu_cxx::__normal_iterator<std::pair<dealii::TriaIterator<dealii::CellAccessor<3,
3> >, unsigned int>* const*,
std::vector<std::pair<dealii::TriaIterator<dealii::CellAccessor<3, 3> >,
unsigned int>*,
std::allocator<std::pair<dealii::TriaIterator<dealii::CellAccessor<3, 3> >,
unsigned int>*> > > > const&>,
boost::_bi::list2<boost::reference_wrapper<dealii::WorkStream::internal::Implementation3::WorkerAndCopier<std::pair<dealii::TriaIterator<dealii::CellAccessor<3,
3> >, unsigned int>*, dealii::internal::DataOut::ParallelData<3, 3>,
dealii::DataOutBase::Patch<3, 3> > >, boost::arg<1> > >,
tbb::auto_partitioner const>::execute()
#18 /user2/hanks318/dealii/dealii_pre/lib/libdeal_II.g.so.8.3.0:
tbb::internal::custom_scheduler<tbb::internal::IntelSchedulerTraits>::local_wait_for_all(tbb::task&,
tbb::task*)
#19 /user2/hanks318/dealii/dealii_pre/lib/libdeal_II.g.so.8.3.0:
tbb::internal::arena::process(tbb::internal::generic_scheduler&)
#20 /user2/hanks318/dealii/dealii_pre/lib/libdeal_II.g.so.8.3.0:
tbb::internal::market::process(rml::job&)
#21 /user2/hanks318/dealii/dealii_pre/lib/libdeal_II.g.so.8.3.0:
tbb::internal::rml::private_worker::run()
#22 /user2/hanks318/dealii/dealii_pre/lib/libdeal_II.g.so.8.3.0:
tbb::internal::rml::private_worker::thread_routine(void*)
--------------------------------------------------------
I can't understand the meaning of distorted image of mapping, non positive
volume etc...
The below picture is the 3D mesh I used for the calculation(it is refined
twice by triangulation.refine_global(2);)
<https://lh3.googleusercontent.com/-mg4mtHj0cPo/V8bL5xBmV3I/AAAAAAAAAAc/AuOBTpAZ5aAnbd3fj4DQWylbPIOh5L2jgCLcB/s1600/mesh_3D.png>
>From this picture, It seems that The Mesh is fine. I don't know what the
problem is.
And FYI, I used "SphericalManifold<dim> boundary" in this code, even though
the boundary is not sphere. I used this to make the boundary smoother after
every refinement.
But, When I didn't use this, The problem is not solved. So, I thought it is
not related to this problem
I also attach the step-20.cc that I revised a little bit.(I commented what
I change in step-20.cc)
please help me to solve this problem....
Thank you
Regards,
Kyusik
--
The deal.II project is located at http://www.dealii.org/
For mailing list/forum options, see
https://groups.google.com/d/forum/dealii?hl=en
---
You received this message because you are subscribed to the Google Groups
"deal.II User Group" group.
To unsubscribe from this group and stop receiving emails from it, send an email
to dealii+unsubscr...@googlegroups.com.
For more options, visit https://groups.google.com/d/optout.
B_coord3D.inp
Description: chemical/gamess-input
/* ---------------------------------------------------------------------
*
* Copyright (C) 2005 - 2015 by the deal.II authors
*
* This file is part of the deal.II library.
*
* The deal.II library is free software; you can use it, redistribute
* it, and/or modify it under the terms of the GNU Lesser General
* Public License as published by the Free Software Foundation; either
* version 2.1 of the License, or (at your option) any later version.
* The full text of the license can be found in the file LICENSE at
* the top level of the deal.II distribution.
*
* ---------------------------------------------------------------------
*
* Author: Wolfgang Bangerth, Texas A&M University, 2005, 2006
*/
// @sect3{Include files}
#include <deal.II/base/quadrature_lib.h>
#include <deal.II/base/logstream.h>
#include <deal.II/base/function.h>
#include <deal.II/lac/block_vector.h>
#include <deal.II/lac/full_matrix.h>
#include <deal.II/lac/block_sparse_matrix.h>
#include <deal.II/lac/solver_cg.h>
#include <deal.II/lac/precondition.h>
#include <deal.II/lac/iterative_inverse.h>
#include <deal.II/grid/tria.h>
#include <deal.II/grid/grid_generator.h>
#include <deal.II/grid/tria_accessor.h>
#include <deal.II/grid/tria_iterator.h>
#include <deal.II/dofs/dof_handler.h>
#include <deal.II/dofs/dof_renumbering.h>
#include <deal.II/dofs/dof_accessor.h>
#include <deal.II/dofs/dof_tools.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/numerics/vector_tools.h>
#include <deal.II/numerics/matrix_tools.h>
#include <deal.II/numerics/data_out.h>
#include <fstream>
#include <iostream>
#include <deal.II/fe/fe_raviart_thomas.h>
#include <deal.II/base/tensor_function.h>
//revised part1~(for input of B_coord3D.inp)
#include <deal.II/grid/grid_in.h>
#include <deal.II/grid/manifold_lib.h>
//revised part1~
namespace Step20
{
using namespace dealii;
template <int dim>
class MixedLaplaceProblem
{
public:
MixedLaplaceProblem (const unsigned int degree);
void run ();
private:
void make_grid_and_dofs ();
void assemble_system ();
void solve ();
//void compute_errors () const; //revised part2 (the exact solution is not exact anymore, so I won't use it)
void output_results () const;
const unsigned int degree;
Triangulation<dim> triangulation;
FESystem<dim> fe;
DoFHandler<dim> dof_handler;
BlockSparsityPattern sparsity_pattern;
BlockSparseMatrix<double> system_matrix;
BlockVector<double> solution;
BlockVector<double> system_rhs;
};
template <int dim>
class RightHandSide : public Function<dim>
{
public:
RightHandSide () : Function<dim>(1) {}
virtual double value (const Point<dim> &p,
const unsigned int component = 0) const;
};
template <int dim>
class PressureBoundaryValues : public Function<dim>
{
public:
PressureBoundaryValues () : Function<dim>(1) {}
virtual double value (const Point<dim> &p,
const unsigned int component = 0) const;
};
/* revised part3 (ExactSolution is not needed any more)
template <int dim>
class ExactSolution : public Function<dim>
{
public:
ExactSolution () : Function<dim>(dim+1) {}
virtual void vector_value (const Point<dim> &p,
Vector<double> &value) const;
};
*/
template <int dim>
double RightHandSide<dim>::value (const Point<dim> &/*p*/,
const unsigned int /*component*/) const
{
return 0;
}
template <int dim>
double PressureBoundaryValues<dim>::value (const Point<dim> &p,
const unsigned int /*component*/) const
{
//const double alpha = 0.3;
//const double beta = 1;
//return -(alpha*p[0]*p[1]*p[1]/2 + beta*p[0] - alpha*p[0]*p[0]*p[0]/6);
return 1+std::abs(p[0]+p[1]); //revised part4 (I just want to make sure x and y can be accessed in this class)
}
/* revised part5(exact sol)
template <int dim>
void
ExactSolution<dim>::vector_value (const Point<dim> &p,
Vector<double> &values) const
{
Assert (values.size() == dim+1,
ExcDimensionMismatch (values.size(), dim+1));
const double alpha = 0.3;
const double beta = 1;
values(0) = alpha*p[1]*p[1]/2 + beta - alpha*p[0]*p[0]/2;
values(1) = alpha*p[0]*p[1];
values(2) = -(alpha*p[0]*p[1]*p[1]/2 + beta*p[0] - alpha*p[0]*p[0]*p[0]/6);
}
*/
template <int dim>
class KInverse : public TensorFunction<2,dim>
{
public:
KInverse () : TensorFunction<2,dim>() {}
virtual void value_list (const std::vector<Point<dim> > &points,
std::vector<Tensor<2,dim> > &values) const;
};
template <int dim>
void
KInverse<dim>::value_list (const std::vector<Point<dim> > &points,
std::vector<Tensor<2,dim> > &values) const
{
Assert (points.size() == values.size(),
ExcDimensionMismatch (points.size(), values.size()));
for (unsigned int p=0; p<points.size(); ++p)
{
values[p].clear ();
for (unsigned int d=0; d<dim; ++d)
//values[p][d][d] = 1.0;
//revised part6 (I just want to check x and y can be accessed in this class)
values[p][d][d] = 1.0+std::abs(points[p][0])+std::abs(points[p][1]);
}
}
template <int dim>
MixedLaplaceProblem<dim>::MixedLaplaceProblem (const unsigned int degree)
:
degree (degree),
fe (FE_RaviartThomas<dim>(degree), 1,
FE_DGQ<dim>(degree), 1),
dof_handler (triangulation)
{}
template <int dim>
void MixedLaplaceProblem<dim>::make_grid_and_dofs ()
{
//revised part7(input 3D Plasma Boundary)
GridIn<dim> grid_in;
grid_in.attach_triangulation(triangulation);
std::ifstream input_file("B_coord3D.inp");
grid_in.read_ucd(input_file);
static const SphericalManifold<dim> boundary;
triangulation.set_all_manifold_ids_on_boundary(0);
triangulation.set_manifold(0,boundary);
//revised part7
triangulation.refine_global (2);
dof_handler.distribute_dofs (fe);
DoFRenumbering::component_wise (dof_handler);
std::vector<types::global_dof_index> dofs_per_component (dim+1);
DoFTools::count_dofs_per_component (dof_handler, dofs_per_component);
const unsigned int n_u = dofs_per_component[0],
n_p = dofs_per_component[dim];
std::cout << "Number of active cells: "
<< triangulation.n_active_cells()
<< std::endl
<< "Total number of cells: "
<< triangulation.n_cells()
<< std::endl
<< "Number of degrees of freedom: "
<< dof_handler.n_dofs()
<< " (" << n_u << '+' << n_p << ')'
<< std::endl;
BlockDynamicSparsityPattern dsp(2, 2);
dsp.block(0, 0).reinit (n_u, n_u);
dsp.block(1, 0).reinit (n_p, n_u);
dsp.block(0, 1).reinit (n_u, n_p);
dsp.block(1, 1).reinit (n_p, n_p);
dsp.collect_sizes ();
DoFTools::make_sparsity_pattern (dof_handler, dsp);
sparsity_pattern.copy_from(dsp);
system_matrix.reinit (sparsity_pattern);
solution.reinit (2);
solution.block(0).reinit (n_u);
solution.block(1).reinit (n_p);
solution.collect_sizes ();
system_rhs.reinit (2);
system_rhs.block(0).reinit (n_u);
system_rhs.block(1).reinit (n_p);
system_rhs.collect_sizes ();
}
template <int dim>
void MixedLaplaceProblem<dim>::assemble_system ()
{
QGauss<dim> quadrature_formula(degree+2);
QGauss<dim-1> face_quadrature_formula(degree+2);
FEValues<dim> fe_values (fe, quadrature_formula,
update_values | update_gradients |
update_quadrature_points | update_JxW_values);
FEFaceValues<dim> fe_face_values (fe, face_quadrature_formula,
update_values | update_normal_vectors |
update_quadrature_points | update_JxW_values);
const unsigned int dofs_per_cell = fe.dofs_per_cell;
const unsigned int n_q_points = quadrature_formula.size();
const unsigned int n_face_q_points = face_quadrature_formula.size();
FullMatrix<double> local_matrix (dofs_per_cell, dofs_per_cell);
Vector<double> local_rhs (dofs_per_cell);
std::vector<types::global_dof_index> local_dof_indices (dofs_per_cell);
const RightHandSide<dim> right_hand_side;
const PressureBoundaryValues<dim> pressure_boundary_values;
const KInverse<dim> k_inverse;
std::vector<double> rhs_values (n_q_points);
std::vector<double> boundary_values (n_face_q_points);
std::vector<Tensor<2,dim> > k_inverse_values (n_q_points);
const FEValuesExtractors::Vector velocities (0);
const FEValuesExtractors::Scalar pressure (dim);
typename DoFHandler<dim>::active_cell_iterator
cell = dof_handler.begin_active(),
endc = dof_handler.end();
for (; cell!=endc; ++cell)
{
fe_values.reinit (cell);
local_matrix = 0;
local_rhs = 0;
right_hand_side.value_list (fe_values.get_quadrature_points(),
rhs_values);
k_inverse.value_list (fe_values.get_quadrature_points(),
k_inverse_values);
for (unsigned int q=0; q<n_q_points; ++q)
for (unsigned int i=0; i<dofs_per_cell; ++i)
{
const Tensor<1,dim> phi_i_u = fe_values[velocities].value (i, q);
const double div_phi_i_u = fe_values[velocities].divergence (i, q);
const double phi_i_p = fe_values[pressure].value (i, q);
for (unsigned int j=0; j<dofs_per_cell; ++j)
{
const Tensor<1,dim> phi_j_u = fe_values[velocities].value (j, q);
const double div_phi_j_u = fe_values[velocities].divergence (j, q);
const double phi_j_p = fe_values[pressure].value (j, q);
local_matrix(i,j) += (phi_i_u * k_inverse_values[q] * phi_j_u
- div_phi_i_u * phi_j_p
- phi_i_p * div_phi_j_u)
* fe_values.JxW(q);
}
local_rhs(i) += -phi_i_p *
rhs_values[q] *
fe_values.JxW(q);
}
for (unsigned int face_no=0;
face_no<GeometryInfo<dim>::faces_per_cell;
++face_no)
if (cell->at_boundary(face_no))
{
fe_face_values.reinit (cell, face_no);
pressure_boundary_values
.value_list (fe_face_values.get_quadrature_points(),
boundary_values);
for (unsigned int q=0; q<n_face_q_points; ++q)
for (unsigned int i=0; i<dofs_per_cell; ++i)
local_rhs(i) += -(fe_face_values[velocities].value (i, q) *
fe_face_values.normal_vector(q) *
boundary_values[q] *
fe_face_values.JxW(q));
}
cell->get_dof_indices (local_dof_indices);
for (unsigned int i=0; i<dofs_per_cell; ++i)
for (unsigned int j=0; j<dofs_per_cell; ++j)
system_matrix.add (local_dof_indices[i],
local_dof_indices[j],
local_matrix(i,j));
for (unsigned int i=0; i<dofs_per_cell; ++i)
system_rhs(local_dof_indices[i]) += local_rhs(i);
}
}
class SchurComplement : public Subscriptor
{
public:
SchurComplement (const BlockSparseMatrix<double> &A,
const IterativeInverse<Vector<double> > &Minv);
void vmult (Vector<double> &dst,
const Vector<double> &src) const;
private:
const SmartPointer<const BlockSparseMatrix<double> > system_matrix;
const SmartPointer<const IterativeInverse<Vector<double> > > m_inverse;
mutable Vector<double> tmp1, tmp2;
};
SchurComplement::SchurComplement (const BlockSparseMatrix<double> &A,
const IterativeInverse<Vector<double> > &Minv)
:
system_matrix (&A),
m_inverse (&Minv),
tmp1 (A.block(0,0).m()),
tmp2 (A.block(0,0).m())
{}
void SchurComplement::vmult (Vector<double> &dst,
const Vector<double> &src) const
{
system_matrix->block(0,1).vmult (tmp1, src);
m_inverse->vmult (tmp2, tmp1);
system_matrix->block(1,0).vmult (dst, tmp2);
}
class ApproximateSchurComplement : public Subscriptor
{
public:
ApproximateSchurComplement (const BlockSparseMatrix<double> &A);
void vmult (Vector<double> &dst,
const Vector<double> &src) const;
void Tvmult (Vector<double> &dst,
const Vector<double> &src) const;
private:
const SmartPointer<const BlockSparseMatrix<double> > system_matrix;
mutable Vector<double> tmp1, tmp2;
};
ApproximateSchurComplement::ApproximateSchurComplement (const BlockSparseMatrix<double> &A)
:
system_matrix (&A),
tmp1 (A.block(0,0).m()),
tmp2 (A.block(0,0).m())
{}
void ApproximateSchurComplement::vmult (Vector<double> &dst,
const Vector<double> &src) const
{
system_matrix->block(0,1).vmult (tmp1, src);
system_matrix->block(0,0).precondition_Jacobi (tmp2, tmp1);
system_matrix->block(1,0).vmult (dst, tmp2);
}
void ApproximateSchurComplement::Tvmult (Vector<double> &dst,
const Vector<double> &src) const
{
vmult (dst, src);
}
template <int dim>
void MixedLaplaceProblem<dim>::solve ()
{
PreconditionIdentity identity;
IterativeInverse<Vector<double> > m_inverse;
m_inverse.initialize(system_matrix.block(0,0), identity);
m_inverse.solver.select("cg");
static ReductionControl inner_control(1500, 0., 1.e-13);
m_inverse.solver.set_control(inner_control);
Vector<double> tmp (solution.block(0).size());
{
Vector<double> schur_rhs (solution.block(1).size());
m_inverse.vmult (tmp, system_rhs.block(0));
system_matrix.block(1,0).vmult (schur_rhs, tmp);
schur_rhs -= system_rhs.block(1);
SchurComplement
schur_complement (system_matrix, m_inverse);
ApproximateSchurComplement
approximate_schur_complement (system_matrix);
IterativeInverse<Vector<double> >
preconditioner;
preconditioner.initialize(approximate_schur_complement, identity);
preconditioner.solver.select("cg");
preconditioner.solver.set_control(inner_control);
SolverControl solver_control (solution.block(1).size(),
1e-12*schur_rhs.l2_norm());
SolverCG<> cg (solver_control);
cg.solve (schur_complement, solution.block(1), schur_rhs,
preconditioner);
std::cout << solver_control.last_step()
<< " CG Schur complement iterations to obtain convergence."
<< std::endl;
}
{
system_matrix.block(0,1).vmult (tmp, solution.block(1));
tmp *= -1;
tmp += system_rhs.block(0);
m_inverse.vmult (solution.block(0), tmp);
}
}
/*
template <int dim>
void MixedLaplaceProblem<dim>::compute_errors () const
{
const ComponentSelectFunction<dim>
pressure_mask (dim, dim+1);
const ComponentSelectFunction<dim>
velocity_mask(std::make_pair(0, dim), dim+1);
ExactSolution<dim> exact_solution;
Vector<double> cellwise_errors (triangulation.n_active_cells());
// As already discussed in step-7, we have to realize that it is
// impossible to integrate the errors exactly. All we can do is
// approximate this integral using quadrature. This actually presents a
// slight twist here: if we naively chose an object of type
// <code>QGauss@<dim@>(degree+1)</code> as one may be inclined to do (this
// is what we used for integrating the linear system), one realizes that
// the error is very small and does not follow the expected convergence
// curves at all. What is happening is that for the mixed finite elements
// used here, the Gauss points happen to be superconvergence points in
// which the pointwise error is much smaller (and converges with higher
// order) than anywhere else. These are therefore not particularly good
// points for integration. To avoid this problem, we simply use a
// trapezoidal rule and iterate it <code>degree+2</code> times in each
// coordinate direction (again as explained in step-7):
QTrapez<1> q_trapez;
QIterated<dim> quadrature (q_trapez, degree+2);
// With this, we can then let the library compute the errors and output
// them to the screen:
VectorTools::integrate_difference (dof_handler, solution, exact_solution,
cellwise_errors, quadrature,
VectorTools::L2_norm,
&pressure_mask);
const double p_l2_error = cellwise_errors.l2_norm();
VectorTools::integrate_difference (dof_handler, solution, exact_solution,
cellwise_errors, quadrature,
VectorTools::L2_norm,
&velocity_mask);
const double u_l2_error = cellwise_errors.l2_norm();
std::cout << "Errors: ||e_p||_L2 = " << p_l2_error
<< ", ||e_u||_L2 = " << u_l2_error
<< std::endl;
}
*/
template <int dim>
void MixedLaplaceProblem<dim>::output_results () const
{
std::vector<std::string> solution_names;
switch (dim)
{
case 2:
solution_names.push_back ("u");
solution_names.push_back ("v");
solution_names.push_back ("p");
break;
case 3:
solution_names.push_back ("u");
solution_names.push_back ("v");
solution_names.push_back ("w");
solution_names.push_back ("p");
break;
default:
Assert (false, ExcNotImplemented());
}
DataOut<dim> data_out;
data_out.attach_dof_handler (dof_handler);
data_out.add_data_vector (solution, solution_names);
data_out.build_patches (degree+1);
std::ofstream output ("solution_new.vtk");
data_out.write_vtk (output);
}
template <int dim>
void MixedLaplaceProblem<dim>::run ()
{
make_grid_and_dofs();
assemble_system ();
solve ();
//compute_errors (); //revised part8 (exact sol part)
output_results ();
}
}
int main ()
{
try
{
using namespace dealii;
using namespace Step20;
deallog.depth_console (0);
MixedLaplaceProblem<3> mixed_laplace_problem(0);
mixed_laplace_problem.run ();
}
catch (std::exception &exc)
{
std::cerr << std::endl << std::endl
<< "----------------------------------------------------"
<< std::endl;
std::cerr << "Exception on processing: " << std::endl
<< exc.what() << std::endl
<< "Aborting!" << std::endl
<< "----------------------------------------------------"
<< std::endl;
return 1;
}
catch (...)
{
std::cerr << std::endl << std::endl
<< "----------------------------------------------------"
<< std::endl;
std::cerr << "Unknown exception!" << std::endl
<< "Aborting!" << std::endl
<< "----------------------------------------------------"
<< std::endl;
return 1;
}
return 0;
}
