Reference documentation for deal.II version 9.5.0
\(\newcommand{\dealvcentcolon}{\mathrel{\mathop{:}}}\) \(\newcommand{\dealcoloneq}{\dealvcentcolon\mathrel{\mkern-1.2mu}=}\) \(\newcommand{\jump}[1]{\left[\!\left[ #1 \right]\!\right]}\) \(\newcommand{\average}[1]{\left\{\!\left\{ #1 \right\}\!\right\}}\)
Loading...
Searching...
No Matches
output.h
Go to the documentation of this file.
1// ---------------------------------------------------------------------
2//
3// Copyright (C) 2010 - 2023 by the deal.II authors
4//
5// This file is part of the deal.II library.
6//
7// The deal.II library is free software; you can use it, redistribute
8// it, and/or modify it under the terms of the GNU Lesser General
9// Public License as published by the Free Software Foundation; either
10// version 2.1 of the License, or (at your option) any later version.
11// The full text of the license can be found in the file LICENSE.md at
12// the top level directory of deal.II.
13//
14// ---------------------------------------------------------------------
15
16
17#ifndef dealii_mesh_worker_output_h
18#define dealii_mesh_worker_output_h
19
20#include <deal.II/base/config.h>
21
25
27
29
30
32
33namespace MeshWorker
34{
35 namespace Assembler
36 {
58 {
59 public:
64
74 void
75 initialize(const unsigned int n_points, const unsigned int n_vectors);
76
81 void
82 initialize_stream(std::ostream &stream);
83
91 template <int dim>
92 void
93 initialize_info(DoFInfo<dim> &info, bool face);
94
98 template <int dim>
99 void
100 assemble(const DoFInfo<dim> &info);
101
105 template <int dim>
106 void
107 assemble(const DoFInfo<dim> &info1, const DoFInfo<dim> &info2);
108
109 private:
114 template <typename T>
115 void
116 write(const T &t) const;
117
123 void
124 write_endl() const;
125
129 unsigned int n_vectors;
133 unsigned int n_points;
134
138 std::ostream *os;
139 };
140
141 //----------------------------------------------------------------------//
142
143 template <typename T>
144 inline void
145 GnuplotPatch::write(const T &d) const
146 {
147 if (os == nullptr)
148 deallog << d;
149 else
150 (*os) << d;
151 }
152
153
154 inline void
156 {
157 if (os == nullptr)
158 deallog << std::endl;
159 else
160 (*os) << std::endl;
161 }
162
163
165 : n_vectors(numbers::invalid_unsigned_int)
166 , n_points(numbers::invalid_unsigned_int)
167 , os(nullptr)
168 {}
169
170
171 inline void
172 GnuplotPatch::initialize(const unsigned int np, const unsigned int nv)
173 {
174 n_vectors = nv;
175 n_points = np;
176 }
177
178
179 inline void
181 {
182 os = &stream;
183 }
184
185
186 template <int dim>
187 inline void
189 {
190 if (face)
191 info.initialize_quadrature(Utilities::fixed_power<dim - 1>(n_points),
192 n_vectors + dim);
193 else
194 info.initialize_quadrature(Utilities::fixed_power<dim>(n_points),
195 n_vectors + dim);
196 }
197
198
199 template <int dim>
200 inline void
202 {
203 const unsigned int np = info.n_quadrature_points();
204 const unsigned int nv = info.n_quadrature_values();
205 const unsigned int patch_dim =
206 (info.face_number == numbers::invalid_unsigned_int) ? dim : (dim - 1);
207 const unsigned int row_length = n_points;
208 // If patches are 1d, end the
209 // patch after a row, else end
210 // it after a square
211 const unsigned int row_length2 =
212 (patch_dim == 1) ? row_length : (row_length * row_length);
213
214 // AssertDimension(np, Utilities::fixed_power<dim>(n_points));
215 AssertDimension(nv, n_vectors + dim);
216
217
218 for (unsigned int k = 0; k < np; ++k)
219 {
220 if (k % row_length == 0)
221 write_endl();
222 if (k % row_length2 == 0)
223 write_endl();
224
225 for (unsigned int i = 0; i < nv; ++i)
226 {
227 write(info.quadrature_value(k, i));
228 write('\t');
229 }
230 write_endl();
231 }
232 }
233
234
235 template <int dim>
236 inline void
238 {
239 assemble(info1);
240 assemble(info2);
241 }
242 } // namespace Assembler
243} // namespace MeshWorker
244
246
247#endif
void write(const T &t) const
Definition output.h:145
void initialize_stream(std::ostream &stream)
Definition output.h:180
void assemble(const DoFInfo< dim > &info)
Definition output.h:201
void initialize_info(DoFInfo< dim > &info, bool face)
Definition output.h:188
void initialize(const unsigned int n_points, const unsigned int n_vectors)
Definition output.h:172
unsigned int face_number
Definition dof_info.h:90
unsigned int n_quadrature_values() const
number & quadrature_value(const unsigned int k, const unsigned int i)
void initialize_quadrature(const unsigned int np, const unsigned int nv)
unsigned int n_quadrature_points() const
#define DEAL_II_NAMESPACE_OPEN
Definition config.h:472
#define DEAL_II_NAMESPACE_CLOSE
Definition config.h:473
#define AssertDimension(dim1, dim2)
LogStream deallog
Definition logstream.cc:37
static const unsigned int invalid_unsigned_int
Definition types.h:213