1
This is ./rheolef.info, produced by makeinfo version 4.13 from
4
This file documents the the `rheolef' system for performing finite
7
Copyright (C) 1996, 1998, 2001, 2002 by Pierre Saramito. All rights
10
Permission is granted to make and distribute verbatim copies of this
11
manual provided the copyright notice and this permission notice are
12
preserved on all copies.
14
Permission is granted to copy and distribute modified versions of this
15
manual under the conditions for verbatim copying, provided that the
16
entire resulting derived work is distributed under the terms of a
17
permission notice identical to this one.
19
INFO-DIR-SECTION Mathematics
21
* rheolef: rheolef The finite element system.
25
File: rheolef.info, Node: Top, Prev: (dir), Up: (dir)
27
Rheolef, the finite element system.
28
***********************************
30
This file documents the the `rheolef' system for performing finite
35
* Abstract:: Preliminary information.
36
* Installing:: How to install `rheolef'.
37
* Problems:: Reporting bugs.
43
* FAQ for developers:: How to develop `rheolef'
44
* Copying:: How you can copy and share `rheolef'.
45
* Concept Index:: Index of concepts.
46
* Program Index:: Index of programs.
47
* Class Index:: Index of classes.
48
* Bilinear Form Index:: Index of bilinear forms.
49
* Approximation Index:: Index of polynomial approximations.
50
* Function Index:: Index of functions.
51
* File Format Index:: Index of file formats.
52
* Related Tool Index:: Index of related tools.
55
File: rheolef.info, Node: Abstract, Up: Top
60
`rheolef' is a computer environment that serves as a convenient
61
_laboratory_ for computations involving finite element-like methods.
62
It provides a set of unix commands and C++ algorithms and containers.
63
This environment is currently under development.
65
Algorithms are related, from one hand to _sparse matrix_ basic linear
66
algebra, e.g. x+y, A*x, A+B, A*B, ... From other hand, we focus
67
development on preconditionned linear and non-linear solver e.g. direct
68
and iterative methods for A*x=b, iterative solvers for Stokes and
69
Bingham flows problems. Future developments will consider also
70
viscoelastic flows problems.
72
Containers covers first the classic graph data structure for sparse
73
matrix formats and finite element _meshes_. An higher level of
74
abstraction is provided by containers related to approximate finite
75
element _spaces_, discrete _fields_ and bilinear _forms_. Current
76
developments concerns distributed sparse matrix. Future developments
77
will consider also 3D anisotropic _adaptative_ mesh generation.
79
Please send all comments and bug reports by electronic mail to
80
<Pierre.Saramito@imag.fr>.
83
File: rheolef.info, Node: Installing, Up: Top
88
(Source file: `configure.ac')
93
Successfull compilation and installation require the GNU `make'
94
command. Check your installation by running a small test:
97
2.1 Reading the documentation
98
=============================
100
`rheolef' comes with three documentations:
103
* an users guide with a set of complete examples
105
* a reference manual for unix commands an C++ classes
107
* the full browsable source code in html format
109
All these documentations are downloadable in various formats from the
110
`rheolef' home page. These files are also locally installed during
111
the `make install' stage. The users guide is generated in `pdf' format:
114
The reference manual is generated in `html' format:
115
firefox http://ljk.imag.fr/membres/Pierre.Saramito/rheolef/rheolef.html
116
and you could add it to your bookmarks. Moreover, after performing
117
`make install', unix manuals are available for commands and classes:
122
The browsable source code is generated in `html' format:
123
firefox http://ljk.imag.fr/membres/Pierre.Saramito/rheolef/source_html/
125
2.2 Using and alternative installation directory
126
================================================
128
The default install prefix is `/usr/local'. If you prefer an
129
alternative directory, e.g. `HOME/sys', you should enter:
130
./configure --prefix=$HOME/sys
131
In that case, you have to add the folowing lines at the end of your
132
`.profile' or `.bash_profile' file:
133
export PATH="$PATH:$HOME/sys/bin"
134
export LD_LIBRARY_PATH="$LD_LIBRARY_PATH:$HOME/sys/lib"
135
export MANPATH="$MANPATH:$HOME/man"
136
assuming you are using the `sh' or the `bash' unix interpreter.
137
Conversely, if you are using `csh' or `tcsh', add at the bottom of
139
set path = ($path $HOME/sys/bin)
140
setenv LD_LIBRARY_PATH "$LD_LIBRARY_PATH:$HOME/sys/lib"
141
setenv MANPATH "$MANPATH:$HOME/sys/man"
142
If you are unure about unix interpreter, get the SHELL value:
144
Then, you have to _source_ the modified file, e.g.:
146
hpux users should replace `LD_LIBRARY_PATH' by `SHLIB_PATH'. If you
147
are unure, get the shared library path variable:
148
rheolef-config --shlibpath-var
149
Check also your installation by compiling a sample program with the
152
Since it is a common mistake, you can later, at each time, test your
153
run-time environment sanity with:
154
rheolef-config --check
156
2.3 Recommanded library
157
=======================
159
The use of implicit numerical schemes leads to the resolution of
160
large sparse linear systems. `rheolef' installation optionally
161
recognizes both `umfpack' and `spooles', two multifrontal direct
162
solvers libraries, and `taucs', an out-of-core sparse direct solver,
163
When no efficient sparse direct solver library is not available,
164
rheolef uses a direct solver based upon a traditionnal skyline Choleski
165
factorisation combined with a reverse Cuthill-Mc Kee reordering.
166
Here is the time complexity for 2d and 3d regular grid-based Poisson
169
factorization N^2 N log N
171
resolution N^3/2 N log N
173
factorization N^7/3 N^2
175
resolution N^5/2 N^4/3
176
The out-of-core strategy is efficient for very large problems or if
177
your computer has few memory. The complexity estimation has been
178
shown by A. George, NESTED DISSECTION OF A REGULAR FINITE ELEMENT
179
MESH, SIAM J. Numer. Anal., 10:345-363, 1973.
181
2.4 The umfpack library
182
=======================
184
The `umfpack' library, if present, is auto-detected and used. It is
185
widely available, as a debian package.
187
2.5 The spooles library
188
=======================
190
The optional `spooles-2.2' multifrontal linear solver from A. George,
191
october 1998, is recommended for solving efficiently large linear
194
It is available at `http://www.netlib.org/linalg/spooles/'. This
195
library is free software. A default skyline strategy is used when
196
it is not founded by `configure'. Both the version 2.2. and the
197
older `2.0' library interface are supported by `rheolef'.
199
First build the spooles library:
202
gzip -dc < ../spooles.2.2.tar.gz | tar xvf -
203
make CC=gcc CFLAGS="-O9" lib
204
Better, you can build spooles as a shared library, as follow
206
make CC=gcc CFLAGS="-O9 -fPIC" lib
210
gcc -shared -o ../libspooles.so *.o
214
Note that on `hpux' systems, the `.sl' suffix is used instead of
215
the `.so' one. Then, go to the rheolef directory:
216
./configure --with-spooles=/usr/local/src/spooles-2.2
220
2.6 The taucs library
221
=====================
223
The optional `taucs-2.0' out-of-core sparse linear solver from S.
224
Toledo and D. Chen, may 2002, is recommended for very large linear
225
systems or if your computer has few memory. `taucs' can factor a
226
matrix whose factors are larger than the main memory by storing the
227
factors on disk files. In that case, this strategy is significantly
228
faster than paging activity when running an in-core algorithm such as
231
It is available at `http://www.tau.ac.il/~stoledo/taucs/'. This
232
library is free software.
234
The `taucs' installation requires `metis', `lapack' and `blas'
235
libraries. Assuming all these libraries are installed in
236
`/usr/local/lib' and the taucs header `taucs.h' is located in
237
`/usr/local/include', the configuration for `rheolef' writes:
239
./configure --with-taucs-ldadd='-ltaucs -llapack -lblas -lmetis -lg2c'
240
An alternative is to `blas' is to use the more efficient `atlas'
241
libraries. See the installation notes for the taucs library.
243
2.7 Too long compilations
244
=========================
246
On platforms with small memory, the `c++' compilation could be long.
247
By default, the `configure' scripts uses a `-O2' compilation level.
248
Turn it off either by:
250
setenv CXXFLAGS "-O0"
253
or by editing directly `config/config.mk'.
255
2.8 Portability issues
256
======================
258
rheolef is based on STL. Therefore you need a compiler that
259
supports the ANSI C++ standards. STL and rheolef especially depend
260
heavily on the template support of the compiler.
262
current version of rheolef has only been tested with
263
GNU C++ compiler (g++-2.95.4, 2.96, 3.0.2, 3.0.4, 4.0, 4.3)
264
Compaq C++ compiler (V6.3, V6.5)
265
on the following operating systems:
266
linux debian, pentium 3 and 4, using GNU C++
267
mac os x (i386-apple-darwin9.8.0), using GNU C++
268
sun solaris 2.8, ultrasparc, using GNU C++
269
compacq OSF1 alpha (V5.1, V4.0), using Compacq C++
270
The GNU C++ compiler is a free software available at
271
`http://www.gnu.org'. The Compacq C++ compiler is available on
272
Compacq (DEC) platforms.
274
Previous versions was tested also on the KAI C++ compiler
275
(KCC-3.2b) and the CRAY C++ compiler (CC-3.0.2.0). on the
276
following operating systems:
277
aix 4.1 ibm sp2, multiprocessor based on rs6000
278
cray unicos t3e, multiprocessor based on 64 bits dec alpha
280
hpux 9.05 and 10.20, hppa
281
linux (redhat, suze and slackware), pentium (2,3,pro)
283
sun solaris 2.5.1, sparc
284
The KAI C++ is available at `http://www.kai.com'. The CRAY C++
285
is available on CRAY platforms. Nevertheless, the current version
286
has no more been tested on these compilers and systems.
288
If you are running a another compiler and operating system combination,
289
please run the non-regression test suite:
291
and reports us the result. If all is ok, please, send us a mail, and
292
we will add your configuration to the list, and else we will try to
293
circumvent the problem.
295
2.9 Known portability limitations
296
=================================
298
On compacq alpha, the old cxx V6.1 (DIGITAL UNIX V4.0) compiles well
299
but scratches on non-regression tests, due to a bug in shared
300
libraries managment. Please, update to cxx V6.3 or more, it will
303
On some hpux and GNU C++ version combinations, building shared
304
libraries could be broken. Uses instead:
305
./configure --disable-shared
306
It should work. Building static libraries leads to larger executable
307
files, and then more disk space available. This limitation depends upon
308
compilers and systems compatibilities, not upon rheolef source code
311
Please read `http://www.gnu.org/software/gcc/install/specific.html'.
312
A recent version of `binutils' and a hpux patch are recommended.
314
Osman F. Buyukisik <<osman.buyukisik@ae.ge.com>> reported that
315
export LIBS=" -ldce -static "
317
works also well on hpux-10.20.
324
Using rheolef suggests the use of the following optional tools:
336
2.11 Advanced configuration
337
===========================
339
Note that each `--with' option has a corresponding `--without'
344
Turns compile-time optimization to maximum available.
347
`--with-bigfloat=DIGITS10'
348
Turn on Bruno Haible's `cln' arbitrary precision float
349
library with DIGITS10 precision,
351
See `http://clisp.cons.org/~haible/packages-cln.html'.
352
Default is off. The DIGITS10 value is optional and
353
default DIGITS10 value is 60. This feature is still
357
Set the home directory for the `cln' library. Default
358
DIR_CLN is `/usr/local/math'.
360
`--with-doubledouble'
361
Uses `doubledouble' class as the default rheolef `Float' type.
362
This option is usefull for a quadruple-like precision on
363
machines where this feature is not or incorrectly
364
available. Default is off.
366
`--enable-debian-packaging=DEBIAN_PACKAGING'
367
Generate files the debian way by respecting the `Debian
370
`--with-boost-incdir=INCDIR_BOOST'
371
Turns on/off the `boost' boost/numeric/uBlas dense and
372
sparse matrix library. This library is required by
376
Produce a `c++ -g' like code.
379
With debugging, also uses dynamic allocation runtime checking
380
with `dmalloc' when available.
382
`--enable-distributed'
383
Turns on/off the distributed memory and parallel computation
384
features. Default is off. This feature is currently in
385
development. When on, configure try to detect either the
386
mpi, and cotch or the parmetis libraries.
389
Turns on/off the distributed memory and parallel computation
390
features. Default is on when distributed is on. This
391
feature is currently in development. When on, configure
392
try to detect either the scotch or the parmetis libraries.
395
`--with-scotch=LIBDIR_SCOTCH'
396
Turns on/off the `scotch' distributed mesh partitionner.
397
Check in LIBDIR_SCOTCH for `libptscotchparmetis.so',
398
`libptscotch.so' and `libptscotcherrexit.so' or
399
corresponding `.a' libraries. Default is to auto-detect
400
when mpi is available.
403
`--with-parmetis=LIBDIR_PARMETIS'
404
Turns on/off the `parmetis' distributed mesh partitionner.
405
Check in LIBDIR_PARMETIS for `libparmetis.so',
406
`libparmetis.a' and also `libmetis.a' `libmetis.so'.
407
Default is to autodetect when mpi is available and scotch
411
`--with-blas-ldadd=RHEO_LDADD_BLAS'
412
Turns on/off the `blas' libraries. Check in
413
RHEO_LDADD_BLAS for `libblas.so', or the corresponding
414
`.a' libraries. Default is to auto-detect. This
415
librariy is required for the pastix library.
418
`--with-pastix=LIBDIR_PASTIX'
419
Turns on/off the `scotch' distributed mesh partitionner.
420
Check in LIBDIR_SCOTCH for `libptscotchparmetis.so',
421
`libptscotch.so' and `libptscotcherrexit.so' or
422
corresponding `.a' libraries. Default is to auto-detect
423
when mpi is available.
425
`--with-umfpack=LIBDIR_UMFPACK'
426
`--with-umfpack-incdir=INCDIR_UMFPACK'
427
Turns on/off the `umfpack' version 4.3 multifrontal direct solver.
428
Default INCDIR_UMFPACK is LIBDIR_UMFPACK. Check in
429
LIBDIR_UMFPACK for `libumfpack.so' or `libumfpack.a' and
430
for header INCDIR_UMFPACK/umfpack.h or
431
INCDIR_UMFPACK/umfpack/umfpack.h. When this library
432
is not available, the direct solver bases upon a
433
traditionnal skyline Choleski factorisation combined with
434
a reverse Cuthill-Mc Kee reordering.
436
`--with-taucs-ldadd=TAUCS_LDADD'
437
`--with-taucs-incdir=TAUCS_INCDIR'
438
Turns on/off the `taucs' version 2.0 out-of-core sparse direct
439
solver. When this library is not available, the
440
configure script check for the spooles multifrontal (see
443
`--with-spooles-libdir=LIBDIR_SPOOLES'
444
`--with-spooles-incdir=INCDIR_SPOOLES'
445
Turns on/off the `spooles' version 2.2 multifrontal direct solver.
446
Default INCDIR_SPOOLES is LIBDIR_SPOOLES. Check in
447
LIBDIR_SPOOLES for `libspooles.so', `libspooles.a' or
448
`spooles.a' and for header INCDIR_SPOOLES/FrontMtx.h.
449
When this library is not available, the direct solver
450
bases upon a traditionnal skyline Choleski factorisation
451
combined with a reverse Cuthill-Mc Kee reordering. Generic or
452
hardware-dependant optimization
454
2.12 Future configure options
455
=============================
457
These options will be implemented in the future:
460
Defines `long double' as the default rheolef `Float' type.
461
Default is off and defines `double' as `Float' type.
462
Warning: most of architectures does not provides
463
mathematical library in `long double', e.g.
464
`sqrt((long double)2)' returns only double precision
465
result on Sun architectures, despites it uses a double memory
466
space. Thus use `--with-doubledouble' instead.
468
2.13 Rebuilding the documentation
469
=================================
471
You can rebuild the documentation by entering:
474
The user's manual goes into `doc/usrman/' and the the reference manual
477
The user's manual becomes available in `.dvi', `.ps' and `.pdf'
480
Regenerating the full documention requires also optional tools.
483
xdvi doc/usrman/usrman.dvi
484
ghostview doc/usrman/usrman.ps
485
xpdf doc/usrman/usrman.pdf
487
The reference manual is available in `.info', `.dvi', `.ps', `.pdf',
488
`.html' and `.txt' formats:
491
info -f doc/refman/rheolef.info
492
xdvi doc/refman/rheolef.dvi
493
ghostview doc/refman/rheolef.ps
494
xpdf doc/refman/rheolef.pdf
495
lynx doc/refman/rheolef_toc.html
496
vi doc/refman/rheolef.txt
498
2.14 Documentation using doxygen
499
================================
501
`doxygen' is a documentation system and can be found at
502
`http://www.doxygen.org'.
504
`doxygen' has built-in support to generate inheritance diagrams for C++
505
classes. `doxygen' can use the `dot' tool from `graphviz' to generate
506
more advanced diagrams and graphs. `graphviz' is an open-sourced,
507
cross-platform graph drawing toolkit from AT&T and Lucent Bell Labs
508
and can be found at `http://www.research.att.com/sw/tools/graphviz/'.
511
mozilla file:`pwd`/doc/doxygen/html/index.html
513
2.15 Compile-time optional tools
514
================================
516
Installation has been tested with the indicated version number.
517
Later version should also work well.
521
The definitive site for this code is
523
`http://www-epidem.plantsci.cam.ac.uk/~kbriggs/doubledouble.html'.
525
The 2.2 version is included in this distribution for convenience.
528
Bruno Haible's arbitrary precision float library.
530
See `http://clisp.cons.org/~haible/packages-cln.html'.
533
For dynamic allocation checking in conjuction with `configure
537
2.16 The directory structure
538
============================
541
portability related files and documentation extractors
544
a collection of useful C++ classes and functions
545
for smart pointers, handling files and
546
directories, using strings, timing.
549
sparse matrix library, renumbering and factorization,
553
finite element library, mesh, spaces and forms.
556
reference manual, user's guide, web site.
559
demonstrations as in users'guide.
561
2.17 Compiling from development version
562
=======================================
564
The installation process is preceeded by an extra `bootstrap' stage,
565
and the `configure' stage takes an extra `--enable-ginac' flag. An
568
./configure --enable-ginac
571
You will need some extra-development tools that are used by
574
* for configuration files:
581
* for automatically builted files:
588
* for maintaining documentation:
595
latex 2e (2001/06/01)
601
* Version numbers are indicative: others versionnumbers sould also work.
602
Nevertheless, these version numbers are known to work well.
604
* For participating to the rheolef concurrent development, please
605
contact also <Pierre.Saramito@imag.fr> for concurrent development.
608
File: rheolef.info, Node: Problems, Up: Top
613
This software is still under development. Please, run the `make check'
614
non-regression tests. Send comments and bug repports to
615
<Pierre.Saramito@imag.fr>. Include the version number, which you can
616
find at the top of this document. Also include in your message the
617
input and the output that the program produced and some comments on the
621
File: rheolef.info, Node: Commands, Up: Top
628
* mkgeo_grid command::
630
* grummp2geo command::
631
* tetgen2geo command::
636
* cemagref2field command::
642
* rheolef-config command::
645
File: rheolef.info, Node: mkgeo_grid command, Up: Commands
647
4.1 mkgeo_grid - build a regular grid mesh in 1d, 2d or 3d
648
==========================================================
650
(Source file: `nfem/bin/mkgeo_grid')
655
mkgeo_grid OPTIONS [NX [NY [NZ]]]
660
The following command build a triangular based 2d 10x10 grid of the
662
mkgeo_grid -t 10 > square-10.geo
664
or in one comand line:
665
mkgeo_grid -t 10 | geo -
670
This command is usefull when testing programs on simple geometries.
671
It avoid the preparation of an input file for a mesh generator. The
672
optional NX, NY and NZ arguments are integer that specifies the
673
subdivision in each direction. By default NX=10, NY=NX and NZ=NY.
674
The mesh files goes on standard output.
676
The command supports all the possible element types: edges, triangles,
677
rectangles, tetraedra, prisms and hexahedra.
686
2d mesh using triangles.
689
2d mesh using quadrangles (rectangles).
692
3d mesh using tetraedra.
695
3d mesh using prisms.
698
3d mesh using hexahedra.
703
The geometry can be any [a,b] segment, [a,b]x[c,d] rectangle or
704
[a,b]x[c,d]x[f,g] parallelotope. By default a=c=f=0 and b=d=g=1, thus,
705
the unit boxes are considered. For instance, the following command
706
meshes the [-2,2]x[-1.5, 1.5] rectangle:
707
mkgeo_grid -t 10 -a -2 -b 2 -c -1.5 -d 1.5 | geo -
720
The boundary domains for boxes uses names: `left', `right',
721
`top', `bottom',`front' and `back'. Instead of splitting the
722
boundary in faces, this option groups all boundary faces in one
723
domain named `boundary'.
725
mkgeo_grid -t 10 -boundary | geo -
731
The whole domain is splitted into two subdomains: `east' and
732
`west', This option is used for testing computations with
733
subdomains (e.g. transmission problem; see the user manual).
735
mkgeo_grid -t 10 -region | geo -
741
The corners (four in 2D and eight in 3D) are defined as OD-domains.
742
This could be usefull for some special boundary conditions.
744
mkgeo_grid -t 10 -corner | geo -
745
mkgeo_grid -T 5 -corner | geo -
747
Coordinate system option
748
------------------------
750
Most of rheolef codes are coordinate-system independant. The
751
coordinate system is specified in the geometry file `.geo'.
753
the 2d mesh is axisymmetric.
756
File: rheolef.info, Node: mfield command, Up: Commands
758
4.2 `mfield' - handle multi-field files
759
=======================================
761
(Source file: `nfem/bin/mfield.cc')
766
mfield {-FIELD-NAME}+ [-IDIR] FILENAME
771
This command performs a _velocity_ graphical output for vector-valued
773
mfield -u square.mfield -velocity
774
mfield -u square.mfield -deformation
775
mfield -s square.mfield -proj -tensor
776
It is also convenient for selecting some scalar fields, e.g. to
777
pipe to another processor
778
mfield -psi square.mfield | field -
779
mfield -u0 square.mfield | field -
780
mfield -u1 square.mfield | field -
785
Read and output multiple finite element fields from file, in field
786
text file format. Fields are preceded by a label, in the following
787
`.mfield' file format:
798
Such labeled file format could be easily generated using the
799
`catchmark' stream manipulator (*note iorheo class::).
800
cout << catchmark("u") << uh
801
<< catchmark("p") << ph
802
<< catchmark("psi") << psih;
808
add DIR to the RHEOPATH search path. See also *note geo
809
class:: for RHEOPATH mechanism.
812
specifies the name of the file containing the input field.
815
all fields are selected for stdout printing.
818
read field on standard input instead on a file.
821
number of digits used to print floating point values when
822
using the `-geo' option. Default depends upon the machine
823
precision associated to the `Float' type.
826
print messages related to graphic files created and command
827
system calls (this is the default).
830
does not print previous messages.
833
Round the output, in order for the non-regression tests to
834
be portable across platforms. Floating point are machine-dependent.
835
The round value is optional. default value is 1e-10.
842
Render vector-valued fields as arrows using `plotmtv',
843
`mayavi' or `vtk'. The the numeric extension is assumed:
844
use e.g. `-u' option instead of `-u0', `-u1' options.
847
Render vector-valued fields as deformed mesh using
848
`plotmtv', `mayavi' or `vtk'. The the numeric extension
849
is assumed: use e.g. `-u' option instead of `-u0', `-u1'
853
Render tensor-valued fields as ellipses using `mayavi'.
854
The the numeric extension is assumed: use e.g. `-tau'
855
option instead of `-tau00', `-tau01'... Warning: in
859
scale vector values by specified factor. Default value is 1.
862
Convert all selected fields to P1-continuous approximation by
863
using a L2 projection. Only valid yet together with the
866
Render tool selection
867
---------------------
870
use `mayavi' tool. This is the default for bi- and
871
tridimensional geometries. The environment variable
872
RHEOLEF_RENDER may be set to one of the choices below to override
876
use `plotmtv' plotting command.
879
use `vtk' tool. Old style for 3d problems: superseted by
883
Others render options
884
---------------------
889
Use (color/gray scale/black and white) rendering. Color
890
rendering is the default.
894
rendering mode: suitable for red-blue anaglyph 3D stereoscopic
898
vector norm isolines intervals are filled with color.
901
vector norm isolines are drawed by using colored mesh edges.
905
print messages related to graphic files created and command
906
system calls (this is the default).
909
does not print previous messages.
912
clear temporary graphic files (this is the default).
915
does not clear temporary graphic files.
918
execute graphic command (this is the default).
921
does not execute graphic command. Generates only graphic
922
files. This is usefull in conjuction with the `-noclean'
928
Check if selected fields appear in file. Otherwise, send some
932
File: rheolef.info, Node: grummp2geo command, Up: Commands
934
4.3 grummp2geo - convert grummp mesh in geo format
935
==================================================
937
(Source file: `nfem/bin/grummp2geo')
942
grummp OPTIONS INPUT[.g] [INPUT[.dmn]]
947
Convert a grummp `.g' into `.geo' one. The output goes to standart
948
output. The `.dmn' file specifies the domain names, since `grummp'
949
mesh generators use numbers as domain labels.
954
grummp2geo toto.g toto.dmn > toto.geo
957
File: rheolef.info, Node: tetgen2geo command, Up: Commands
959
4.4 tetgen2geo - convert tetgen mesh in geo format
960
==================================================
962
(Source file: `nfem/bin/tetgen2geo')
967
tetgen2geo OPTIONS INPUT[.node] INPUT[.ele] INPUT[.face] INPUT[.dmn]
972
Convert a tetgen mesh into `.geo' one. The output goes to standart
973
output. The `.dmn' file specifies the domain names, since `tetgen'
974
mesh generator uses numbers as domain labels.
980
tetgen2geo toto.1.node toto.1.ele toto.1.face toto.dmn > toto.geo
985
This auxilliary `.dmn' file defines the boundary domain names as used
986
by Rheolef, since `tetgen' uses numeric labels for domains.
993
Default is to output a version 2 `.geo' file format. *Note geo
994
command::. With the `-noupgrade', a version 1 file format
998
File: rheolef.info, Node: geo command, Up: Commands
1000
4.5 `geo' - plot a finite element mesh
1001
======================================
1003
(Source file: `nfem/bin/xgeo.cc')
1008
geo OPTIONS MESH[.geo]
1013
Plot a 2d or 3d finite element mesh.
1015
The input file format is described with the `geo' class manual (*note
1023
geo -vtk -fill bielle.geo
1025
geo -plotmtv square.geo
1027
Input file specification
1028
------------------------
1031
specifies the name of the file containing the input mesh.
1032
The ".geo" extension is assumed.
1035
add GEODIR to the mesh search path. This mechanism
1036
initializes a search path given by the environment variable
1037
`RHEOPATH'. If the environment variable `RHEOPATH' is not
1038
set, the default value is the current directory
1039
(*note geo class::).
1042
read mesh on standard input instead on a file.
1045
when mesh comes from standard input, the mesh name is not
1046
known and is set to "output" by default. This option
1047
allows to change this default. Useful when dealing with
1048
output formats (graphic, format conversion) that creates
1049
auxilliary files, based on this name.
1055
print the number of elements in the mesh and exit.
1058
print the number of nodes in the mesh and exit.
1061
print the smallest edge length in the mesh and exit.
1064
print the largest edge length in the mesh and exit.
1067
print the lower-left-bottom vertice in the mesh and exit.
1070
print the upper-right-top vertice in the mesh and exit.
1075
The default render could also specified by using the RHEOLEF_RENDER
1076
environment variable.
1078
use Mayavi. This is the default.
1087
use Visualization Toolkit.
1090
use x3d visualization tool.
1093
use PlotM atom representation (from lassp, Cornell).
1099
use colors. This is the default.
1102
Option only available with `mayavi'.
1105
Rendering mode suitable for red-blue anaglyph 3D stereoscopic
1106
glasses. Option only available with `mayavi'.
1109
rendering mode. This is the default.
1112
rendering mode. Fill mesh faces using light effects when available.
1115
rendering mode. Shrink 3d elements (with vtk only).
1118
rendering mode. All edges appears as tubes (with vtk only).
1121
rendering mode. Vertices appears as balls (with vtk only).
1124
rendering mode. All edges appears (with vtk only).
1127
rendering mode. cut by plane and clip (with vtk only).
1128
Works with `-full', `-tube' or `-ball' modes. Does not
1129
clip mesh with `-fill' and `-grid' modes. Does not work
1130
yet with `-shrink' mode. This option is still in
1134
set the origin point of the cutting plane when using with
1135
the `-cut' option. Default is (0.5,0.5,0.5).
1137
`-normal NX [NY [NZ]]'
1138
set the normal vector of the cutting plane. when using
1139
with the `-cut' option. Default is (1,0,0).
1142
split each edge by inserting a node. The resulting mesh
1145
Output format options
1146
---------------------
1149
output mesh on standard output stream in geo text file format,
1150
instead of plotting it.
1153
output mesh on standard output stream in bamg text file format,
1154
instead of plotting it.
1157
output mesh on standard output stream in gmsh ascii file format,
1158
instead of plotting it.
1161
output mesh on stream in tegen text file format, instead
1162
of plotting it (TODO).
1165
output mesh on standard output stream in mmg3d text file format,
1166
instead of plotting it.
1169
output mesh on standard output stream in cemagref text file format,
1170
instead of plotting it. Notices that the mesh does not
1171
contains the topographic elevation. If you want to include
1172
the elevation, in the output, use the `field' command with
1173
the `-cemagref' option.
1176
edges and faces are numbered. This numbering is required
1177
for the quadratic approximation.
1180
number of digits used to print floating point values when
1181
using the `-geo' option. Default depends upon the machine
1182
precision of the `Float' type.
1185
Avoid variation of output float formats when using high
1186
precision `doubledouble' or `bigfloat' classes.
1187
This option is used for non-regression test purpose.
1190
Input format options
1191
--------------------
1194
read the mesh in the `.geo' format. This is the default.
1197
read the mesh in the `.bamg' format. Since edges are not
1198
numbered in `.bamg' format, this file format is not
1199
suitable for using P2 elements. This format is supported
1200
only for format conversion purpose. See `-upgrade' for
1201
edge numbering while converting to `.geo'.
1204
read the mesh in the `.msh' format. Since edges in 2D and
1205
faces in 3D are not numbered in `.msh' format, this file
1206
format is not suitable for using P2 elements. This format
1207
is supported only for format conversion purpose. See
1208
`-upgrade' for edge numbering while converting to `.geo'.
1211
read the mesh as the concatenation of `.node', `.ele' and `.face'
1212
tetgen format. Since in 3D are not all numbered in tetgen
1213
format, this file format is not suitable for using P2
1214
elements. This format is supported only for format
1215
conversion purpose. See `-upgrade' for edge numbering
1216
while converting to `.geo'.
1219
read the mesh in the `.mmg3d' format. Since edges and
1220
faces are not numbered in `.mmg3d' format, this file
1221
format is not suitable for using P2 elements. This format
1222
is supported only for format conversion purpose. See
1223
`-upgrade' for edge numbering while converting to `.geo'.
1226
read the mesh on standard input stream in grummp text file format.
1227
The `.g' grummp file format is defined grummp
1228
`.template' specification file
1231
"2" nverts ncells nbfaces
1235
Conversely, for tridimensionnal meshes, the grummp
1236
`.template' specification file is:
1239
"3" nverts ncells nbfaces
1243
Such files can be generated with the `tri' and `tetra'
1244
grummp mesh file generators. The grummp mesh generators suport
1245
multi-regions geometries. These regions are translated to
1246
bi- or tri-dimensional domains, that are no boundary
1247
domains. A file format conversion writes:
1248
geo -input-grummp -geo -upgrade -noverbose - < hex-multi.g > hex-multi.geo
1251
output mesh on standard output stream in qmg text file format.
1252
instead of plotting it. Since edges are not numbered
1253
in `.qmg' format, this file format is not suitable for
1254
using P2 elements. This format is supported only for
1255
format conversion purpose. See `-upgrade' for edge
1256
numbering while converting to `.geo'.
1258
`-dmn DOMAIN-NAMES-FILE'
1259
An auxilliary `.dmn' file defining the boundary domain names
1260
as used by `rheolef' could be used, since `bamg', `mmg3d',
1261
`tetgen' and `grummp' use numeric labels for
1262
domains. Example of a `multi-domain.dmn' with two regions
1263
and four boundary domains file:
1274
`geo' recognized this format as an extension at the end
1275
of the `.bamg', `.mesh' (mmg3d) or the `.g' file. The
1276
complete file format conversion writes:
1277
cat multi-domain.g multi-domain.dmn | \\
1278
geo -input-grummp -upgrade -geo - > multi-domain.geo
1279
If this file is not provided, boundary domains are named
1280
`boundary1', `boundary2', ... and region domain are
1281
named `region1', `region2', ... The
1282
tridimensionnal domain name file contains the
1283
`FaceDomainNames' and `VolumeDomainNames' keywords, for
1284
boundary and region domains, respectively.
1287
read the mesh in the `.cemagref' format. Notices that the
1288
reader skip the topographic elevation informations. If
1289
you want to load the elevation, use the `field' command
1290
with the `-input-cemagref' option.
1297
print messages related to graphic files created and
1298
command system calls (this is the default).
1301
does not print previous messages.
1304
clear temporary graphic files (this is the default).
1307
does not clear temporary graphic files.
1310
execute graphic command (this is the default).
1313
does not execute graphic command. Generates only graphic
1314
files. This is usefull in conjuction with the "-noclean"
1319
used for debug purpose.
1324
You can generate a `postscript', `latex' or `Fig' output of
1325
your mesh, using the `geo' command with `-gnuplot -noclean'
1326
options, and then modifying the `.plot' generated file.
1331
This is the default mesh file format. It contains two entitiess, namely
1332
a mesh and a list of domains. The mesh entity starts with the header.
1333
The header contains the `mesh' keyword followed by a line containing a
1334
format version number, the space dimension (e.g. 1, 2 or 3), the number
1335
of vertices, the number of elements. For file format version 2, when
1336
dimension is three, the number of faces is specified, and then, when
1337
dimension is two or three, the number of edges is also specified.
1338
Follows the vertex coordinates, the elements and the faces or edges.
1339
One element starts with a letter indicating the element type
1360
Then, we have the vertex indexes. A vertex index is numbered in the
1361
C-style, i.e. the first index is numbered 0. An edge is a couple of
1364
Geo version 1 file format (obsolete) does neither give any information
1365
about edges for 2d meshes nor for faces in 3d. This information is
1366
requiered for maintaining P2 approximation file in a consistent way,
1367
since edge and face numbering algorithm should vary from one version to
1368
another, and thus may be stored on file.
1370
Nevertheless, geo version 1 file format is still supported, since it is
1371
a convenient way to build and import meshes (*note geo command::).
1373
A sample mesh is in version 1 writes
1383
the same in version 2
1398
The second entity is a list of domains, that finishes with the end of
1399
file. A domain starts with the `domain' keyword, followed by a domain
1400
name. Then, the format version number (i.e. 1 today), and the number
1401
of sides. Follows the sides
1413
File: rheolef.info, Node: bamg2geo command, Up: Commands
1415
4.6 bamg2geo - convert bamg mesh in geo format
1416
==============================================
1418
(Source file: `nfem/bin/bamg2geo')
1423
bamg2geo OPTIONS INPUT[.bamg] INPUT[.dmn]
1424
bamg2geo OPTIONS INPUT[.bamg] -Cl DOMLABEL
1425
bamg2geo OPTIONS INPUT[.bamg] {-dom DOMNAME}*
1430
Convert a bamg `.bamg' into `.geo' one. The output goes to standart
1431
output. The `.dmn' file specifies the domain names, since `bamg' mesh
1432
generator uses numbers as domain labels.
1437
bamg -g toto.bamgcad -o toto.bamg
1438
bamg2geo toto.bamg toto.dmn > toto.geo
1443
This file describe the boundary of the mesh geometry. A basic example
1444
writes (See bamg documentation for more);
1445
MeshVersionFormatted
1467
This auxilliary `.dmn' file defines the boundary domain names as used
1468
by Rheolef, since `bamg' uses numeric labels for domains.
1476
The domain name file can also specify additional vertices domain
1477
----------------------------------------------------------------
1491
Vertice domain names are usefull for some special boundary conditions.
1498
Default is to output a version 2 `.geo' file format. *Note geo
1499
command::. With the `-noupgrade', a version 1 file format
1502
`-dom DOM1 ... -dom DOMN'
1503
`-Cl {poi|sed|sec|tro|con|NONE}'
1504
Predefined domain name convention. See next section.
1507
File: rheolef.info, Node: mesh2geo command, Up: Commands
1509
4.7 mesh2geo - convert mmg3d mesh in geo format
1510
===============================================
1512
(Source file: `nfem/bin/mesh2geo')
1517
mesh2geo OPTIONS INPUT[.mesh] INPUT[.dmn]
1522
Convert a mesh `.mesh' into `.geo' one. The output goes to standart
1523
output. The `.dmn' file specifies the domain names, since `mmg3d' mesh
1524
generator uses numbers as domain labels.
1529
mmg3d -O 1 -in toto.mesh -out toto-opt.mesh
1530
mesh2geo toto-opt.mesh toto.dmn > toto-opt.geo
1535
This auxilliary `.dmn' file defines the boundary domain names as used
1536
by Rheolef, since `mmg3d' uses numeric labels for domains.
1543
Default is to output a version 2 `.geo' file format. *Note geo
1544
command::. With the `-noupgrade', a version 1 file format
1548
File: rheolef.info, Node: field command, Up: Commands
1550
4.8 `field' - plot a field
1551
==========================
1553
(Source file: `nfem/bin/xfield.cc')
1558
field [-IDIR] FILENAME
1563
Read and output a finite element field from file.
1569
field square.field -black-and-white
1571
field box.field -volume
1573
Input file specification
1574
------------------------
1577
add DIR to the RHEOPATH search path. See also *note geo
1578
class:: for RHEOPATH mechanism.
1581
specifies the name of the file containing the input field.
1584
read field on standard input instead on a file.
1587
when field data comes from standard input, the field name
1588
is not known and is set to "output" by default. This
1589
option allows to change this default. Useful when dealing
1590
with output formats (graphic, format conversion) that
1591
creates auxilliary files, based on this name.
1594
Input format options
1595
--------------------
1598
read the data in the `.field' format. This is the default.
1601
input data in cemagref file format. The field represents
1602
the topographic elevation and is assumed to be a `P1'
1605
Output format options
1606
---------------------
1610
output field on standard output stream in field text file format.
1613
output field on standard output stream in bamg file format.
1616
output field on standard output stream in mmg3d file format.
1619
output field on standard output stream in gmsh file format.
1622
output field on standard output stream in cemagref file format.
1623
The field represents the topographic elevation and is
1624
assumed to be a `P1' approximation. This option is
1625
suitable for file format conversion:
1628
number of digits used to print floating point values.
1629
Default depends upon the machine precision associated to the
1633
Round the output, in order for the non-regression tests to
1634
be portable across platforms. Floating point are machine-dependent.
1637
File format conversions
1638
-----------------------
1640
Input and output options can be combined for efficient file format
1642
field -input-cemagref -field - < mytopo.cemagref > mytopo.field
1643
field myfile.field -bamg > myfile.bb
1644
can be performed in one pass:
1645
field -input-cemagref -bamg - < mytopo.cemagref > myfile.bb
1652
print the min (resp. max) value of the scalar field and then exit.
1657
The default render could also specified by using the RHEOLEF_RENDER
1658
environment variable.
1660
use `gnuplot' tool. This is the default for
1661
monodimensional geometries.
1664
use `mayavi' tool. This is the default for bi- and
1665
tridimensional geometries.
1668
use `plotmtv' plotting command, for bidimensional geometries.
1671
use `vtk' tool. This is the old (no more supported) tool
1672
for tridimensional geometries. Could be used in 2d also.
1681
Use (color/gray scale/black and white) rendering. Color
1682
rendering is the default.
1686
rendering mode: suitable for red-blue anaglyph 3D stereoscopic
1690
isoline intervals are filled with color. This is the
1694
draw isolines by using lines.
1697
mesh edges appear also.
1700
prevent mesh edges drawing. This is the default.
1703
Convert all selected fields to P1-continuous approximation by
1704
using a L2 projection. Usefull when input is a
1705
discontinuous field (e.g. P0 or P1d).
1708
extract an iso-value set of points (1d), lines (2d) or
1709
surfaces (3d) for a specified value. Output the surface in
1710
text format (with `-text' option) or using a graphic
1711
driver. This option requires the `vtk' code.
1714
do not draw isosurface.
1718
This is an experimental volume representation by using ray cast
1719
(mayavi and vtk graphics).
1722
For 2D visualizations, the isovalue table contains
1723
regularly spaced values from fmin to fmax, the bounds of
1726
`-n-iso-negative INT'
1727
The isovalue table is splitted into negatives and positives values.
1728
Assume there is n_iso=15 isolines: if 4 is requested by
1729
this option, then, there will be 4 negatives isolines,
1730
regularly spaced from fmin to 0 and 11=15-4 positive
1731
isolines, regularly spaced from 0 to fmax. This
1732
option is usefull when plotting e.g. vorticity or stream
1733
functions, where the sign of the field is representative.
1737
Write or do not write labels for values on isolines on 2d
1738
contour plots when using `plotmtv' output.
1741
For two dimensional field, represent values as elevation
1742
in the third dimension.
1745
Prevent from the `elevation' representation. This is the
1749
applies a multiplicative factor to the field. This is
1750
useful e.g. in conjonction with the `elevation' option.
1751
The default value is 1.
1754
show a cut by a specified plane. The cutting plane is
1755
specified by its origin point and normal vector. This
1756
option requires the `vtk' code.
1758
`-origin FLOAT [FLOAT [FLOAT]]'
1759
set the origin of the cutting plane. Default is (0.5,
1762
`-normal FLOAT [FLOAT [FLOAT]]'
1763
set the normal of the cutting plane. Default is (1, 0, 0).
1770
print messages related to graphic files created and command
1771
system calls (this is the default).
1774
does not print previous messages.
1777
clear temporary graphic files (this is the default).
1780
does not clear temporary graphic files.
1783
execute graphic command (this is the default).
1786
does not execute graphic command. Generates only graphic
1787
files. This is usefull in conjuction with the `-noclean'
1793
It contains a header and a list values at degrees of freedom. The
1794
header contains the `field' keyword followed by a line containing a
1795
format version number (presently 1), the number of degrees of freedom
1796
(i.e. the number of values listed), the mesh file name without the
1797
`.geo' extension the approximation (e.g. P1, P2, etc), and finaly the
1798
list of values: A sample field file (compatible with the sample mesh
1799
example presented in command manual; see *note geo command::) writes:
1812
field cube.field -cut -normal 0 1 0 -origin 0.5 0.5 0.5 -vtk
1813
This command send to `vtk' the cutted 2d plane of the 3d field.
1814
field cube.field -cut -normal 0 1 0 -origin 0.5 0.5 0.5 -text > cube-cut.field
1815
This command generates the cutted 2d field and its associated mesh.
1816
field cube.field -iso 0.5 -plotmtv
1817
This command draws the isosurface.
1818
field cube.field -iso 0.5 -text > isosurf.geo
1819
This command generates the isosurface as a 3d surface mesh in
1820
`.geo' format. This is suitable for others treatments.
1823
File: rheolef.info, Node: cemagref2field command, Up: Commands
1825
4.9 cemagref2field - convert cemagref topography in field and geo formats
1826
=========================================================================
1828
(Source file: `nfem/bin/cemagref2field')
1833
cemagref2field OPTIONS INPUT[.cemagref]
1838
Convert a cemagref `.cemagref' topography data file into `.field' and
1839
its associated `.geo' files. The `.field' output goes to standart
1840
output, while the compressed `.geo' is created as INPUT.`geo.gz' using
1846
cemagref2field toto.cemagref > toto-z.field
1847
cemagref2field -cemagref toto.dat > toto-z.field
1848
and a `toto.geo.gz'.
1854
specifies directly the input file, for an alternate suffix
1858
File: rheolef.info, Node: qmg2geo command, Up: Commands
1860
4.10 qmg2geo - convert qmg mesh in geo format
1861
=============================================
1863
(Source file: `nfem/bin/qmg2geo')
1868
Convert a qmg `.qmg' into `.geo' one. The output goes to standart
1874
qmg2geo [-qmg] INPUT.qmg INPUT.dmn
1879
This auxilliary `.dmn' file defines the boundary domain names as used
1880
by Rheolef, since `qmg' uses numeric labels for domains.
1889
File: rheolef.info, Node: space command, Up: Commands
1891
4.11 `space' - print a finite element space
1892
===========================================
1894
(Source file: `nfem/bin/xspace.cc')
1899
`space' FILENAME[.space]
1904
Read and re-output a finite element space.
1915
This is the default space file format. The file starts with
1916
the header: the "space" keyword followed by a line containing a
1917
format version number (i.e. 1 today), and the number of degrees of
1918
freedom. Follows the numbering and an optional 'B' label, if the
1919
degree of freedom is blocked. A sample space is:
1930
File: rheolef.info, Node: msh2geo command, Up: Commands
1932
4.12 msh2geo - convert gmsh mesh in geo format
1933
==============================================
1935
(Source file: `nfem/bin/msh2geo')
1940
msh2geo OPTIONS INPUT[.msh]
1945
Convert a gmsh `.msh' into `.geo' one. The output goes to standart
1951
gmsh -2 toto.mshcad -o toto.msh
1952
msh2geo toto.msh > toto.geo
1953
See the `gmsh' documentation for a detailed description of the
1954
`.mshcad' input file for `gmsh'.
1957
File: rheolef.info, Node: branch command, Up: Commands
1959
4.13 `branch' - handle a family of fields
1960
=========================================
1962
(Source file: `nfem/bin/branch.cc')
1967
branch [OPTIONS] FILENAME
1972
Generates vtk file colection for visualization with paraview:
1973
branch output.branch -paraview
1978
Read and output a branch of finite element fields from file, in field
1981
Input file specification
1982
------------------------
1985
add DIR to the RHEOPATH search path. See also *note geo
1986
class:: for RHEOPATH mechanism.
1989
specifies the name of the file containing the input field.
1992
read field on standard input instead on a file.
1995
Number of digits used to print floating point values when
1996
using the `-geo' option. Default depends upon the machine
1997
precision associated to the `Float' type.
2000
Output and render specification
2001
-------------------------------
2004
Extract the i-th record in the file. The output is a field
2005
or multi-field file format.
2008
Output on stdout in `.branch' format. This is the default.
2011
Generate a collection of `vtk' files for using `paraview'.
2014
Generate a single `vtk' file with numbered fields.
2017
Run 1d animation using `gnuplot'.
2020
This driver is unsupported for animations.
2026
`-topography FILENAME[.field[.gz]]'
2027
performs a tridimensionnal elevation view based on the
2031
performs a `P1' projection on the fly. This option is
2032
useful when rendering `P0' data while `vtk' render
2033
requieres `P1' description.
2036
For two dimensional field, represent values as elevation in the
2037
third dimension. This is the default.
2040
Prevent from the elevation representation.
2043
applies a multiplicative factor to the field. This is
2044
useful e.g. in conjonction with the `elevation' option.
2045
The default value is 1.
2048
print messages related to graphic files created and command
2049
system calls (this is the default).
2052
does not print previous messages.
2055
clear temporary graphic files (this is the default).
2058
does not clear temporary graphic files.
2061
execute graphic command (this is the default).
2064
does not execute graphic command. Generates only graphic
2065
files. This is usefull in conjuction with the `-noclean'
2072
The `.branch' file format bases on the `.field' one:
2073
EXAMPLE GENERAL FORM
2077
1 1 11 <version> <nfield=1> <nvalue=N>
2078
time u <key> <field name>
2080
#time 3.14 #<key> <key value 1>
2086
#time 6.28 #<key> <key value N>
2090
The key name is here `time', but could be any string (without
2091
spaces). The previous example contains one `field' at each time
2092
step. Labels appears all along the file to facilitate direct jumps
2093
and field and step skips.
2095
The format supports several fields, such as (t,u(t),p(t)), where u could
2096
be a multi-component (e.g. a vector) field:
2118
File: rheolef.info, Node: cad command, Up: Commands
2120
4.14 `cad' - plot a boundary file
2121
=================================
2123
(Source file: `nfem/bin/cad.cc')
2128
Manipulates geometry boundary files (CAD, Computer Aid Design):
2129
visualizations and file conversions. This command is under
2135
col < FILE[.qmgcad] | cad -input-qmg - OPTIONS...
2141
col < ../data/logo2.qmgcad | cad -input-qmg - -noclean
2144
col < ../data/test2-out.qmgcad | cad -input-qmg -
2145
In the `geomview' window, enter '2as' key sequence to switch to a
2151
The `qmg' demonstration files comes from PC-DOS and contains line-feed
2152
characters that are not filtered by `flex' when reading file. Thus,
2153
the `col' filter is required.
2155
Included surfaces, such as cracks, are hiden in 3D with `geomview'.
2156
Thus, `vtk' rendering with some transparency factor would be better.
2158
`bamg' and `grummp' boundary files are not yet supported. They will
2161
Input file specification
2162
------------------------
2165
specifies the name of the file containing the input mesh.
2166
The ".cad" extension is assumed.
2169
add CADDIR to the mesh search path. This mechanism
2170
initializes a search path given by the environment variable
2171
`RHEOPATH'. If the environment variable `RHEOPATH' is not
2172
set, the default value is the current directory
2173
(*note cad class::).
2176
read mesh on standard input instead on a file.
2178
Input format options
2179
--------------------
2182
read the boundary in the `.cad' format. This is the
2183
default (TODO: not yet, uses qmg instead).
2186
read the boundary in the `.qmgcad' format. In `qmg', this
2187
file format is known as _brep_, that stands for boundary
2190
Output format options
2191
---------------------
2194
output boundary on standard output stream in cad text file format,
2195
instead of plotting it.
2198
use geomview rendering tool. This is the default for
2199
tridimensional data.
2202
use gnuplot rendering tool. This is the default for
2209
Subdivide each Bezier patch in DEGREE*NSUB linear
2210
sub-element. For rendering tools that does not support
2211
Bezier patches, such as `gnuplot' or `geomview'.
2212
Default is NSUB=3 for tridimensional data and NSUB=50 for
2216
Subdivide each Bezier patch in NSUB instead of DEGREE*NSUB.
2219
Subdivide each Bezier patch in of DEGREE*NSUB. This is
2226
print messages related to graphic files created and
2227
command system calls (this is the default).
2230
does not print previous messages.
2233
clear temporary graphic files (this is the default).
2236
does not clear temporary graphic files.
2239
execute graphic command (this is the default).
2242
does not execute graphic command. Generates only graphic
2243
files. This is usefull in conjuction with the "-noclean"
2248
File: rheolef.info, Node: rheolef-config command, Up: Commands
2250
4.15 rheolef-config - get installation directories
2251
==================================================
2253
(Source file: `rheolef-config.in')
2258
The following command returns the rheolef libraries directory:
2259
rheolef-config --libdir
2260
An environment sanity check writes:
2261
rheolef-config --check
2266
This command is usefull when linking executables with rheolef:
2267
libraries locations are required by the link editor. Such
2268
directories are defined while configuring rheolef, before to compile
2269
and install *note Installing::. The `rheolef-config' command returns
2272
Note that `rheolef-config' could be used in Makefiles for the
2273
determination of linker flags.
2275
Another usefull feature is the `--check' option. When `rheolef' is
2276
installed in a user directory, i.e. not as root, the sane run-time
2277
environment depends upon two environment variables. The first one is
2278
the PATH: `bkindir' directory may be present in PATH. The second
2279
environment variable is related to shared libraries, and its name is
2280
system-dependent, e.g. LD_LIBRARY_PATH on most platforms and
2281
SHLIB_PATH on HP-UX. Its content may contains `bindir'.
2282
rheolef-config --shlibpath-var
2283
Since it is a common mistake to have incorrect values for these
2284
variable, for novice users or for adanced ones, especialy when
2285
dealing with several installed versions, the environment sanity check
2287
rheolef-config --check
2288
If there is mistakes, a hint is suggested to fix it and the return
2289
status is 1. Instead, the return status is 0.
2298
print option summary and exit.
2301
install architecture-independent files location.
2304
architecture-dependent files location.
2307
include header directory.
2310
executables directory.
2313
man documentation directory.
2316
object code libraries directory.
2320
read-only architecture-independent data location.
2323
read-only architecture-independent data location; specific for
2327
include compiler flags.
2330
library compiler flags.
2333
the shared library path variable.
2335
`--library-interface-version'
2336
the library interface version.
2338
`--hardcode-libdir-flag-spec'
2339
flag to hardcode a libdir into a binary during linking.
2342
File: rheolef.info, Node: Classes, Up: Top
2351
* form_diag_manip class::
2354
* form_manip class::
2363
* permutation class::
2371
* rheostream class::
2374
File: rheolef.info, Node: cad class, Up: Classes
2376
5.1 `cad' - the boundary definition class
2377
=========================================
2379
(Source file: `nfem/lib/cad.h')
2384
The `cad' class defines a container for the description of the
2385
boundary of a geometry, by using Bezier patches.
2387
The aim is to handle high-level polynomial interpolation together
2388
with curved boundaries (non-polynomial). So, the description bases on
2389
Bezier curves and surfaces.
2391
Also, the adaptive mesh loop requieres interpolation on the boundary,
2392
and this class should facilitates such procedures.
2394
This class is actually under development.
2401
File format conversion
2402
----------------------
2404
Under development. Conversion to geomview and vtk for visualization.
2405
Also to and from qmg, grummp, modulef/ghs, opencascade for
2411
class cad : public smart_pointer<cad_rep> {
2415
typedef cad_rep::size_type size_type;
2421
cad (const std::string& file_name);
2425
size_type dimension() const;
2426
size_type n_face() const;
2428
const point& xmin() const;
2429
const point& xmax() const;
2431
void eval (const cad_element& S, const point& ref_coord, point & real_coord) const;
2432
point eval (const cad_element& S, const point& ref_coord) const;
2436
friend std::istream& operator >> (std::istream& is, cad& x);
2437
friend std::ostream& operator << (std::ostream& os, const cad& x);
2441
File: rheolef.info, Node: geo class, Up: Classes
2443
5.2 `geo' - the finite element mesh class
2444
=========================================
2446
(Source file: `nfem/lib/geo.h')
2451
The `geo' class defines a container for a finite element mesh. This
2452
describes the nodes coordinates and the connectivity. `geo' can
2453
contains domains, usefull for boundary condition setting.
2458
A sample usage of the geo class writes
2461
cout << mayavi << full << omega;
2466
The empty constructor makes an empty geo. A string argument, as
2468
will recursively look for a `square.geo[.gz]' file in the directory
2469
mentionned by the RHEOPATH environment variable, while `gzip'
2470
decompression is assumed. If the file starts with `.' as `./square' or
2471
with a `/' as in `/home/oscar/square', no search occurs. Also, if the
2472
environment variable RHEOPATH is not set, the default value is the
2475
Input and output on streams are available, and manipulators works for
2476
text or graphic output (*note geo command::).
2478
Finally, an STL-like interface provides efficient accesses to nodes,
2479
elements and domains.
2484
The `geo_adapt' functions performs a mesh adaptation to improve the
2485
`P1'-Lagrange interpolation of the `field' given in argument (*note
2488
Axisymetric geometry
2489
--------------------
2491
The `coordinate_system' and `set_coordinate_system' members supports
2492
both `cartesian', `rz' and `zr' (axisymmetric) coordinate systems. This
2493
information is used by the `form' class (*note form class::).
2495
Access to connectivity
2496
----------------------
2498
The folowing code prints triangle vertex numbers
2499
geo omega ("circle");
2500
for (geo::const_iterator i = g.begin(); i != i.end(); i++) {
2501
const geo_element& K = *i;
2502
if (K.name() != 't') continue;
2503
for (geo::size_type j = 0; j < 3; j++)
2504
cout << K [j] << " ";
2507
*Note geo_element internal::.
2509
Access to vertice coordinates
2510
-----------------------------
2512
The folowing code prints vertices coordinate
2513
for (geo::const_iterator_vertex i = g.begin_vertex(); i != g.end_node(); i++) {
2514
const point& xi = *i;
2515
for (geo::size_type j = 0; j < g.dimension(); j++)
2516
cout << xi [j] << " ";
2523
The following code prints edges on domain:
2524
for (geo::const_iterator_domain i = g.begin_domain(); i != i.end_domain(); i++) {
2525
const domain& dom = *i;
2526
if (dom.dimension() != 2) continue;
2527
for (domain::const_iterator j = dom.begin(); j < dom.end(); j++) {
2528
const geo_element& E = *j;
2529
cout << E [0] << " " << E[1] << endl;
2533
*Note domain class::.
2538
RHEOPATH: search path for geo data file. Also *note form class::.
2543
class geo : public smart_pointer<georep> {
2548
typedef georep::plot_options plot_options;
2549
void write_gnuplot_postscript_options (std::ostream& plot, const plot_options& opt) const;
2551
typedef georep::iterator iterator;
2552
typedef georep::const_iterator const_iterator;
2553
typedef georep::elemlist_type elemlist_type;
2554
typedef georep::nodelist_type nodelist_type;
2555
typedef georep::size_type size_type;
2556
typedef georep::domlist_type domlist_type;
2558
typedef georep::const_iterator_node const_iterator_node;
2559
typedef georep::const_iterator_vertex const_iterator_vertex;
2560
typedef georep::const_iterator_domain const_iterator_domain;
2561
typedef georep::iterator_node iterator_node;
2562
typedef georep::iterator_vertex iterator_vertex;
2563
typedef georep::iterator_domain iterator_domain;
2565
// allocators/deallocators:
2567
explicit geo (const std::string& filename, const std::string& coord_sys = "cartesian");
2568
geo(const geo& g, const domain& d);
2569
geo(const nodelist_type& p, const geo& g);
2572
friend geo geo_adapt (const class field& criteria, const Float& hcoef,
2573
bool reinterpolate_criteria = false);
2575
friend geo geo_adapt (const class field& criteria,
2576
const adapt_option_type& = adapt_option_type(),
2577
bool reinterpolate_criteria = false);
2579
friend geo geo_metric_adapt (const field& mh,
2580
const adapt_option_type& = adapt_option_type());
2584
friend std::istream& operator >> (std::istream&, geo&);
2585
friend std::ostream& operator << (std::ostream&, const geo&);
2587
void use_double_precision_in_saving();
2588
std::ostream& dump (std::ostream& s = std::cerr) const;
2589
std::ostream& put_mtv_domains (std::ostream&, size_type=0) const;
2593
const point& vertex (size_type i) const;
2594
const geo_element& element (size_type K_idx) const;
2595
Float measure (const geo_element& K) const;
2597
std::string name() const;
2598
// Family name plus number
2599
std::string basename() const;
2600
// For moving boundary problems
2601
std::string familyname() const;
2602
// Number of moving boundary mesh
2603
size_type number() const;
2604
// Refinment iteration for the current mesh number
2605
size_type version() const;
2606
size_type dimension() const;
2607
size_type map_dimension() const;
2608
std::string coordinate_system () const; // "cartesian", "rz", "zr"
2609
fem_helper::coordinate_type coordinate_system_type() const;
2611
size_type serial_number() const;
2613
const domain& boundary() const;
2614
void build_subregion(const domain& start_from, const domain& dont_cross,
2615
std::string name, std::string name_of_complement="");
2616
// Builds a new domain on the side of domain `dont_cross' on which `start_from'
2617
// lies. These must have no intersection.
2619
size_type size() const;
2620
size_type n_element() const; // same as size()
2621
size_type n_vertex() const;
2622
size_type n_vertice() const;
2623
size_type n_node() const; // same as n_vertex()
2624
size_type n_edge() const ;
2625
size_type n_face() const ;
2626
size_type n_triangle() const ;
2627
size_type n_quadrangle() const ;
2628
size_type n_volume() const ;
2629
size_type n_tetraedra() const ;
2630
size_type n_prism() const ;
2631
size_type n_hexaedra() const ;
2632
size_type n_subgeo(size_type d) const ;
2633
size_type n_element(reference_element::enum_type t) const;
2637
const point& xmin() const;
2638
const point& xmax() const;
2640
meshpoint hatter (const point& x, size_type K) const;
2641
point dehatter (const meshpoint& S) const;
2642
point dehatter (const point& x_hat, size_type e) const;
2644
const_iterator begin() const;
2645
const_iterator end() const;
2646
const_iterator_node begin_node() const;
2647
const_iterator_node end_node() const;
2648
const_iterator_vertex begin_vertex() const;
2649
const_iterator_vertex end_vertex() const;
2652
bool localize (const point& x, geo_element::size_type& element) const;
2653
void localize_nearest (const point& x, point& y, geo_element::size_type& element) const;
2654
bool trace (const point& x0, const point& v, point& x, Float& t, size_type& element) const;
2656
// access to domains:
2657
const domain& get_domain(size_type i) const;
2658
size_type n_domain() const;
2659
bool has_domain (const std::string& domname) const;
2660
const domain& operator[] (const std::string& domname) const;
2661
const_iterator_domain begin_domain() const;
2662
const_iterator_domain end_domain() const;
2664
point normal (const geo_element& S) const;
2665
point normal (const geo_element& K, georep::size_type side_idx) const;
2667
sort_interface(const domain&, const interface&) const;
2669
normal (const class space& Nh, const domain& d,
2670
const std::string& region="") const;
2672
normal (const class space& Nh, const domain& d, const interface& bc) const;
2674
tangent (const class space& Nh, const domain& d,
2675
const std::string& region="") const;
2677
tangent (const class space& Nh, const domain& d, const interface& bc) const;
2678
// Gives normal to d in discontinuous space Nh (P0 or P1d) and can
2679
// initialize orientation of domain through `bc' data structure.
2680
// Currently limited to straight-sided elements.
2682
tangent_spline (const space& Nh, const domain& d, const interface& bc) const;
2684
plane_curvature (const space& Nh, const domain& d,
2685
const interface& bc) const;
2686
// Gives curvature of d in the (x,y) or (r,z) plane.
2687
// Currently limited to straight-sided elements.
2689
plane_curvature_spline (const space& Nh, const domain& d,
2690
const interface& bc) const;
2691
// Gives curvature of d in the (x,y) or (r,z) plane based on a spline interpolation.
2693
plane_curvature_quadratic (const space& Nh, const domain& d,
2694
const interface& bc) const;
2695
// Gives curvature of d in the (x,y) or (r,z) plane based on a local quadratic interpolation.
2698
axisymmetric_curvature (const class space& Nh, const domain& d) const;
2700
axisymmetric_curvature (const class space& Nh, const domain& d,
2701
const interface& bc) const;
2702
// Specific to "rz" and "zr" coordinate systems:
2703
// Gives curvature of d in a plane orthogonal to the (r,z) plane.
2704
// Can also initialize orientation of domain through `bc' data structure.
2705
// Currently limited to straight-sided elements.
2707
interface_process (const domain& d, const interface& bc,
2708
geometric_event& event) const;
2709
// Detects event along domain d sorted according to bc's.
2711
// construction of jump-interface domain data:
2712
void jump_interface(const domain& interface, const domain& subgeo,
2713
std::map<size_type, tiny_vector<size_type> >& special_elements,
2714
std::map<size_type,size_type>& node_global_to_interface) const;
2718
bool operator == (const geo&) const;
2719
bool operator != (const geo&) const;
2723
// build from a list of vertex and elements:
2724
template <class ElementContainer, class VertexContainer>
2725
void build (const ElementContainer&, const VertexContainer&);
2727
void set_name (const std::string&);
2728
void set_coordinate_system (const std::string&); // "cartesian", "rz", "zr"
2730
void label_interface(const domain& dom1, const domain& dom2, const std::string& name);
2733
void erase_domain (const std::string& name);
2734
void insert_domain (const domain&);
2736
// accessors to internals:
2737
bool localizer_initialized () const;
2738
void init_localizer (const domain& boundary, Float tolerance=-1, int list_size=0) const;
2739
const size_type* get_geo_counter() const;
2740
const size_type* get_element_counter() const;
2743
int gnuplot2d (const std::string& basename,
2744
plot_options& opt) const;
2746
// check consistency
2752
friend class spacerep;
2755
void may_have_version_2() const; // fatal if not !
2757
// modifiers to internals:
2759
void set_dimension (size_type d);
2762
iterator_node begin_node();
2763
iterator_node end_node();
2764
iterator_vertex begin_vertex();
2765
iterator_vertex end_vertex();
2769
File: rheolef.info, Node: form_diag_manip class, Up: Classes
2771
5.3 `form_diag_manip' - concatenation on a product space
2772
========================================================
2774
(Source file: `nfem/lib/form_diag_manip.h')
2779
build a form_diag on a product space.
2784
The following example builds a 3-blocks diagonal form:
2791
I_manip << size(3) << Iv << Iv << Iv;
2797
class form_diag_manip {
2802
typedef form_diag::size_type size_type;
2804
// constructors/destructors:
2810
size_type nbloc() const;
2811
size_type nform() const;
2812
form_diag& bloc(size_type i);
2813
const form_diag& bloc(size_type i) const;
2817
form_diag_manip& operator << (const form_diag& a);
2818
form_diag_manip& operator >> (form_diag& a);
2819
friend form_diag_manip& verbose (form_diag_manip &fm);
2820
friend form_diag_manip& noverbose (form_diag_manip &fm);
2821
form_diag_manip& operator<< (form_diag_manip& (*pf)(form_diag_manip&));
2826
friend form_diag_manip& operator << (form_diag_manip &fm, const class sized& s);
2832
std::vector<form_diag> _a;
2836
File: rheolef.info, Node: space class, Up: Classes
2838
5.4 `space' - piecewise polynomial finite element space
2839
=======================================================
2841
(Source file: `nfem/lib/space.h')
2846
The `space' class contains some numbering for unknowns and blocked
2847
_degrees of freedoms_ related to a given mesh and polynomial
2850
Degree of freedom numbering
2851
---------------------------
2853
Numbering of degrees of freedom (or shortly, "dof") depends upon the
2854
mesh and the piecewise polynomial approximation used. This numbering
2855
is then used by the `field' and `form' classes. See also *note field
2856
class:: and *note form class::.
2858
The degree-of-freedom (or shortly, "dof") follows some simple rules,
2859
that depends of course of the approximation. These rules are suitable
2860
for easy `.field' file retrieval (*note field class::).
2864
dof numbers follow element numbers in mesh.
2867
dof numbers follow vertice numbers in mesh.
2870
dof numbers related to vertices follow vertice numbers in mesh.
2871
Follow dof numbers related to edges, in the edge numbering
2875
dof numbers follow element numbers in mesh. In each element, we
2876
loop on vertices following the local vertice order provided in
2879
Unknown and blocked degree of freedom
2880
-------------------------------------
2882
A second numbering is related to the "unknown" and "blocked" degrees of
2884
geo omega("square");
2885
space V(omega,"P1");
2886
V.block ("boundary");
2887
Here, all degrees of freedom located on the domain `"boundary"' are
2893
File output format for `space' is not of practical interest. This file
2894
format is provided for completness. Here is a small example
2905
where the `B' denotes blocked degres of freedom. The method
2906
`set_dof(K, dof_array)' gives the numbering in `dof_array', supplying
2907
an element `K'. The method `index(dof)' gives the unknown or blocked
2908
numbering from the initial numbering, suppling a given `dof'
2913
The method `xdof(K, i)' gives the geometrical location of the i-th
2914
degree of freedom in the element `K'.
2916
Product spaces, for non-scalar fields, such as vector-valued or
2917
tensor-valued (symmetric) fields,
2922
space U (omega, "P1", "vector");
2923
space T (omega, "P1", "tensor");
2924
Then, the number of component depends upon the geometrical dimension of
2927
Arbitrarily product spaces can also be constructed
2929
Spaces for fields defined on a boundary domain can also be constructed
2930
space W (omega, omega["top"], "P1");
2931
The correspondance between the degree-of-freedom numbering for the
2932
space `W' and the global degree-of-freedom numbering for the space `V'
2933
is supported. The method `set_global_dof(S, dof_array)' gives the
2934
global numbering in `dof_array', supplying a boundary element `S'. The
2935
method `domain_index(global_dof)' gives the unknown or blocked
2936
numbering from the initial numbering, suppling a given `global_dof',
2937
related to the space `V'.
2942
class space : public smart_pointer<spacerep> {
2946
typedef spacerep::size_type size_type;
2948
// allocator/deallocator:
2951
space (const geo& g, const std::string& approx_name, const std::string& valued = "scalar");
2952
space (const geo& g, const std::string& approx_name,
2953
const domain& interface, const domain& subgeo, const std::string& valued = "scalar");
2954
// For spaces with a discontinuity through domain `interface', but otherwise continuous
2955
space (const geo& g, const domain& d, const std::string& approx_name);
2956
space(const const_space_component&);
2958
space operator = (const const_space_component&);
2959
friend space operator * (const space&, const space&);
2960
friend space pow (const space&, size_type);
2964
void block () const;
2965
void block_Lagrange () const;
2966
void block (size_type i_comp) const;
2967
void block_Lagrange (size_type i_comp) const;
2968
void block (const std::string& d_name) const;
2969
void block (const domain& d) const;
2970
void block (const std::string& d_name, size_type i_comp) const;
2971
void block (const domain& d, size_type i_comp) const;
2972
void block_dof (size_type dof_idx, size_type i_comp) const;
2973
void block_dof (size_type dof_idx) const;
2974
void unblock_dof (size_type dof_idx, size_type i_comp) const;
2975
void unblock_dof (size_type dof_idx) const;
2976
void set_periodic (const domain& d1, const domain& d2) const;
2977
void set_periodic (const std::string& d1_name, const std::string& d2_name) const;
2978
void lock_components (const domain& d, point locked_dir) const;
2979
void lock_components (const std::string& d_name, point locked_dir) const;
2980
template <class Function>
2981
void lock_components (const std::string& d_name, Function locked_dir) const;
2982
// Allows to enforce vector-boundary conditions.
2983
/* Allows to enforce vector-boundary conditions.
2985
*! Current limitation: only 2D vector fields with components having same FE approximation.
2986
*! The space has a dof for the direction normal to locked_dir, while the value in the
2987
*! direction locked_dir is blocked.
2988
*! The field::at function implements the reconstruction of the original cartesian components.
2990
template <class Function>
2991
void lock_components (const domain& d, Function locked_dir) const;
2995
size_type size() const;
2996
size_type degree () const;
2997
size_type n_blocked() const;
2998
size_type n_unknown() const;
2999
size_type dimension() const;
3000
std::string coordinate_system() const;
3001
fem_helper::coordinate_type coordinate_system_type() const;
3004
const geo& get_geo () const;
3005
const basis& get_basis (size_type i_component = 0) const;
3006
const numbering& get_numbering (size_type i_component = 0) const;
3007
std::string get_approx(size_type i_component = 0) const;
3008
const basis& get_p1_transformation () const;
3010
std::string get_valued() const;
3011
fem_helper::valued_field_type get_valued_type() const;
3013
size_type n_component() const;
3014
space_component operator [] (size_type i_comp);
3015
const_space_component operator [] (size_type i_comp) const;
3017
size_type size_component (size_type i_component = 0) const;
3018
size_type n_unknown_component (size_type i_component = 0) const;
3019
size_type n_blocked_component (size_type i_component = 0) const;
3020
size_type start (size_type i_component = 0) const;
3021
size_type start_unknown (size_type i_component = 0) const;
3022
size_type start_blocked (size_type i_component = 0) const;
3024
size_type index (size_type degree_of_freedom) const;
3025
size_type period_association (size_type degree_of_freedom) const;
3026
bool is_blocked (size_type degree_of_freedom) const;
3027
bool is_periodic (size_type degree_of_freedom) const;
3029
bool has_locked () const;
3030
// for (2D) vector spaces only, field is only blocked along locked_dir
3031
bool is_locked (size_type degree_of_freedom) const;
3032
// for tangential-normal representation
3033
size_type locked_with (size_type degree_of_freedom) const;
3034
size_type index_locked_with (size_type degree_of_freedom) const;
3035
// dof with the other component (thru space::index() as well)
3036
Float locked_component (size_type dof, size_type i_comp) const;
3037
// i-th component of locked direction
3038
Float unlocked_component (size_type dof, size_type i_comp) const;
3039
// i-th component of unlocked direction
3041
void freeze () const;
3042
bool is_frozen() const;
3045
// informations related to boundary spaces
3046
bool is_on_boundary_domain() const;
3047
const domain& get_boundary_domain() const;
3048
const geo& get_global_geo () const;
3049
size_type global_size () const;
3050
size_type domain_dof (size_type global_dof) const;
3051
size_type domain_index (size_type global_dof) const;
3053
// get indices of degrees of freedom associated to an element
3054
// the tiny_vector is a "local index" -> "global index" table
3055
// "global" variants differ only for boundary spaces (non-global uses domain-indices, global use
3056
// mesh-wide indices)
3057
void set_dof (const geo_element& K, tiny_vector<size_type>& idx) const;
3058
void set_global_dof (const geo_element& K, tiny_vector<size_type>& idx) const;
3059
void set_dof (const geo_element& K, tiny_vector<size_type>& idx, size_type i_comp) const;
3060
void set_global_dof (const geo_element& K, tiny_vector<size_type>& idx, size_type i_comp) const;
3061
void set_dof (const geo_element& K, tiny_vector<size_type>& idx, size_type i_comp, reference_element::dof_family_type family) const;
3062
point x_dof (const geo_element& K, size_type i_local) const;
3063
// Hermite/Argyris: returns 1 if the global dof contains the derivative along outward normal of K, -1 else
3064
void dof_orientation(const geo_element& K, tiny_vector<int>& orientation) const;
3066
// piola transformation
3067
meshpoint hatter (const point& x, size_type K) const;
3068
point dehatter (const meshpoint& S) const;
3069
point dehatter (const point& x_hat, size_type e) const;
3073
bool operator == (const space&) const;
3074
bool operator != (const space&) const;
3078
friend std::ostream& operator << (std::ostream&, const space&);
3079
friend std::istream& operator >> (std::istream&, space&);
3080
std::ostream& dump(std::ostream& s = std::cerr) const;
3085
static size_type inquire_size(const geo&, const numbering&);
3089
space(smart_pointer<spacerep>);
3090
space (const space&, const space&); // X*Y
3091
space (const space&, size_type i_comp); // X[i_comp]
3096
File: rheolef.info, Node: form_diag class, Up: Classes
3098
5.5 `form_diag' - diagional form class
3099
======================================
3101
(Source file: `nfem/lib/form_diag.h')
3106
Implement a form using two diagonal matrix.
3111
Here is the computation of the mass matrix for the P1 finite element:
3114
form_diag m(V,"mass");
3126
typedef basic_diag<Float>::size_type size_type;
3128
// allocator/deallocator:
3131
form_diag (const space& X);
3132
form_diag (const space& X, const char* op_name);
3133
form_diag (const class field& dh);
3137
const space& get_space() const;
3138
const geo& get_geo() const;
3139
size_type size() const;
3140
size_type nrow() const;
3141
size_type ncol() const;
3143
// linear algebra (partial):
3145
friend form_diag operator * (const form_diag&, const form_diag&);
3146
friend class field operator * (const form_diag& a, const class field& x);
3147
friend form_diag operator / (const Float& lambda, const form_diag& m);
3149
friend form_diag operator * (const Float& lambda, const form_diag&);
3150
friend form_diag operator + (const form_diag&, const form_diag&);
3151
friend form_diag operator - (const form_diag&);
3152
friend form_diag operator - (const form_diag&, const form_diag&);
3157
friend std::ostream& operator << (std::ostream& s, const form_diag& a);
3163
basic_diag<Float> uu;
3164
basic_diag<Float> bb;
3168
File: rheolef.info, Node: form_manip class, Up: Classes
3170
5.6 `form_manip' - form concatenation on a product space
3171
========================================================
3173
(Source file: `nfem/lib/form_manip.h')
3178
build a form on a product space
3183
a general implementation without switch and with a loop.
3192
typedef form::size_type size_type;
3194
// constructors/destructor:
3200
size_type nform() const;
3201
size_type nrowbloc() const;
3202
size_type ncolbloc() const;
3203
form& bloc(size_type i, size_type j);
3204
const form& bloc(size_type i, size_type j) const;
3208
form_manip& operator << (const form&);
3209
form_manip& operator >> (form&);
3210
form_manip& operator << (form_manip& (*)(form_manip&));
3212
friend form_manip& verbose (form_manip&);
3213
friend form_manip& noverbose (form_manip&);
3218
friend form_manip& operator << (form_manip&, const class size&);
3223
size_type _nrowbloc;
3224
size_type _ncolbloc;
3225
std::vector<form> _a;
3229
File: rheolef.info, Node: field class, Up: Classes
3231
5.7 `field' - piecewise polynomial finite element field
3232
=======================================================
3234
(Source file: `nfem/lib/field.h')
3239
Store degrees of freedom associated to a mesh and a piecewise
3240
polynomial approximation, with respect to the numbering defined by
3241
the underlying *note space class::.
3243
This class contains two vectors, namely unknown and blocked degrees
3244
of freedoms, and the associated finite element space. Blocked and
3245
unknown degrees of freedom can be set by using domain name indexation:
3246
geo omega_h ("circle");
3247
space Vh (omega_h, "P1");
3248
Vh.block ("boundary");
3250
uh ["boundary"] = 0;
3251
Interpolation of a function `u' in a field `uh' with respect to the
3252
interpolation writes:
3253
Float u (const point&);
3254
uh = interpolate (Vh, u);
3259
Here is a complete example using the field class:
3260
Float u (const point& x) { return x[0]*x[1]; }
3262
geo omega_h ("square");
3263
space Vh (omega_h, "P1");
3265
uh = interpolate (Vh, u);
3266
cout << plotmtv << u;
3272
Algebra, such as x+y, x*y, x/y, lambda*x, ...are supported
3274
Transformations applies to all values of a field:
3275
field vh = compose(fabs, uh);
3276
field wh = compose(atan2, uh, vh);
3277
The composition supports also general unary and binary
3280
Vector-valued and tensor-valued field support is yet partial only.
3281
This feature will be more documented in the future.
3284
File: rheolef.info, Node: trace class, Up: Classes
3286
5.8 `trace' - restriction to boundary values operator
3287
=====================================================
3289
(Source file: `nfem/lib/trace.h')
3294
The `trace' operator restricts a field to its values located to a
3295
boundary value domain. The `trace' class is a derived class of the
3296
`form' class (*note form class::). Thus, all operations, such as
3297
linear algebra, that applies to `form', applies also to `trace'.
3302
The following example shows how to plot the trace of a field, from
3305
geo omega ("square");
3306
space Vh (omega "P1");
3307
space Wh (omega, "P1", omega ["top"]);
3308
trace gamma (Vh, Wh);
3320
class trace : public form
3324
// allocator/deallocator:
3327
trace(const space& V, const space& W);
3331
const space& get_space() const;
3332
const space& get_boundary_space() const;
3333
const geo& get_geo() const;
3334
const geo& get_boundary_geo() const;
3338
File: rheolef.info, Node: domain class, Up: Classes
3340
5.9 `domain' - part of a finite element mesh
3341
============================================
3343
(Source file: `nfem/lib/domain.h')
3348
The `domain' class defines a container for a part of a finite element
3349
mesh. This describes the connectivity of edges or faces. This
3350
class is usefull for boundary condition setting.
3355
class domain : public Vector<geo_element> {
3360
typedef Vector<geo_element> sidelist_type;
3361
typedef Vector<geo_element>::size_type size_type;
3363
// allocators/deallocators:
3365
domain(size_type sz = 0, const std::string& name = std::string());
3366
domain(const domain&);
3370
const std::string& name() const;
3371
size_type dimension() const;
3375
void set_name(const std::string&);
3376
void set_dimension(size_type d);
3377
domain& operator=(const domain&);
3378
domain& operator += (const domain&);
3379
friend domain operator + (const domain&, const domain&);
3380
void resize(size_type n);
3381
void cat(const domain& d);
3382
template <class IteratorPair>
3383
void set(IteratorPair p, size_type n, const std::string& name);
3384
template <class IteratorElement>
3385
void set(IteratorElement p, size_type n, size_type dim, const std::string& name);
3389
friend std::ostream& operator << (std::ostream& s, const domain& d);
3390
friend std::istream& operator >> (std::istream& s, domain& d);
3391
std::ostream& put_vtk (std::ostream& vtk, Vector<point>::const_iterator first_p,
3392
Vector<point>::const_iterator last_p) const;
3393
std::ostream& dump (std::ostream& s = std::cerr) const;
3404
File: rheolef.info, Node: form class, Up: Classes
3406
5.10 `form' - representation of a finite element operator
3407
=========================================================
3409
(Source file: `nfem/lib/form.h')
3414
The form class groups four sparse matrix, associated to a bilinear
3415
form on finite element space. This class is represented by four
3416
sparse matrix, associated to unknown and blocked degrees of freedom
3417
of origin and destination spaces (*note space class::).
3422
The operator A associated to a bilinear form a(.,.) by the relation
3423
(Au,v) = a(u,v) could be applied by using a sample matrix notation
3424
A*u, as shown by the following code:
3427
form a(V,V,"grad_grad");
3431
cout << plotmtv << y;
3433
This block-matrix operation is equivalent to:
3434
y.u = a.uu*x.u + a.ub*x.b;
3435
y.b = a.bu*x.u + a.bb*x.b;
3444
typedef csr<Float>::size_type size_type;
3446
// allocator/deallocator:
3449
form (const space& X, const space& Y);
3450
// locked_boundaries means that special vector BCs are applied,
3451
// see space::lock_components()
3452
form (const space& X, const space& Y, const std::string& op_name,
3453
bool locked_boundaries=false);
3454
form (const space& X, const space& Y, const std::string& op_name, const domain& d);
3455
// Currently weighted forms in P2 spaces use a quadrature formula which is
3456
// limited to double precision floats.
3457
form (const space& X, const space& Y, const std::string& op_name, const field& wh,
3458
bool use_coordinate_system_weight=true);
3459
form (const space& X, const space& Y, const std::string& op_name,
3460
const domain& d, const field& wh, bool use_coordinate_system_weight=true);
3461
form (const space& X, const std::string& op_name);
3462
form (const form_diag& X);
3466
const geo& get_geo() const;
3467
const space& get_first_space() const;
3468
const space& get_second_space() const;
3469
size_type nrow() const;
3470
size_type ncol() const;
3471
bool for_locked_boundaries() const { return _for_locked_boundaries; }
3475
Float operator () (const class field&, const class field&) const;
3476
friend class field operator * (const form& a, const class field& x);
3477
field trans_mult (const field& x) const;
3478
friend form trans(const form&);
3479
friend form operator * (const Float& lambda, const form&);
3480
friend form operator + (const form&, const form&);
3481
friend form operator - (const form&);
3482
friend form operator - (const form&, const form&);
3483
friend form operator * (const form&, const form&);
3484
friend form operator * (const form&, const form_diag&);
3485
friend form operator * (const form_diag&, const form&);
3489
friend std::ostream& operator << (std::ostream& s, const form& a);
3490
friend class form_manip;
3496
bool _for_locked_boundaries;
3503
form form_nul(const space& X, const space& Y) ;
3506
File: rheolef.info, Node: branch class, Up: Classes
3508
5.11 `branch' - (t,uh
3509
=====================
3511
(Source file: `nfem/lib/branch.h')
3516
Stores a `field' together with its associated parameter value: a
3517
`branch' variable represents a pair (t,uh(t)) for a specific value of
3518
the parameter t. Applications concern time-dependent problems and
3519
continuation methods.
3521
This class is convenient for file inputs/outputs and building
3522
graphical animations.
3532
This class is under development.
3534
The `branch' class store pointers on `field' class without reference
3535
counting. Thus, `branch' automatic variables cannot be returned by
3536
functions. `branch' variable are limited to local variables.
3538
At each step, meshes are reloaded and spaces are rebuilted, even when
3541
The `vtk' render does not support mesh change during the loop (e.g.
3542
when using adaptive mesh process).
3547
class branch : public std::vector<std::pair<std::string,field> > {
3553
branch (const std::string& parameter_name, const std::string& u_name);
3554
branch (const std::string& parameter_name, const std::string& u_name, const std::string& p_name);
3559
const Float& parameter () const;
3560
const std::string& parameter_name () const;
3561
size_type n_value () const;
3562
size_type n_field () const;
3566
void set_parameter (const Float& value);
3570
// get/set current value
3571
friend std::istream& operator >> (std::istream&, branch&);
3572
friend std::ostream& operator << (std::ostream&, const branch&);
3574
__branch_header header ();
3575
__const_branch_header header () const;
3576
__const_branch_finalize finalize () const;
3578
__obranch operator() (const Float& t, const field& u);
3579
__obranch operator() (const Float& t, const field& u, const field& p);
3581
__iobranch operator() (Float& t, field& u);
3582
__iobranch operator() (Float& t, field& u, field& p);
3585
File: rheolef.info, Node: tensor class, Up: Classes
3587
5.12 `tensor' - a N*N tensor, N=1,2,3
3588
=====================================
3590
(Source file: `nfem/basis/tensor.h')
3595
The `tensor' class defines a 3*3 tensor, as the value of a tensorial
3596
valued field. Basic algebra with scalars, vectors of R^3 (i.e. the
3597
`point' class) and `tensor' objects are supported.
3605
typedef size_t size_type;
3606
typedef Float element_type;
3610
tensor (const Float& init_val = 0);
3611
tensor (Float x[3][3]);
3612
tensor (const tensor& a);
3616
tensor& operator = (const tensor& a);
3617
tensor& operator = (const Float& val);
3621
void fill (const Float& init_val);
3623
void set_row (const point& r, size_t i, size_t d = 3);
3624
void set_column (const point& c, size_t j, size_t d = 3);
3628
Float& operator()(size_type i, size_type j);
3629
Float operator()(size_type i, size_type j) const;
3630
point row(size_type i) const;
3631
point col(size_type i) const;
3632
size_t nrow() const; // = 3, for template matrix compatibility
3633
size_t ncol() const;
3637
friend std::istream& operator >> (std::istream& in, tensor& a);
3638
friend std::ostream& operator << (std::ostream& out, const tensor& a);
3639
std::ostream& put (std::ostream& s, size_type d = 3) const;
3643
friend bool operator == (const tensor&, const tensor&);
3644
friend tensor operator - (const tensor&);
3645
friend tensor operator + (const tensor&, const tensor&);
3646
friend tensor operator - (const tensor&, const tensor&);
3647
friend tensor operator * (Float, const tensor&);
3648
friend tensor operator * (const tensor& a, Float k);
3649
friend tensor operator / (const tensor& a, Float k);
3650
friend point operator * (const tensor& a, const point& x);
3651
friend point operator * (const point& yt, const tensor& a);
3652
point trans_mult (const point& x) const;
3653
friend tensor trans (const tensor& a, size_type d = 3);
3654
friend tensor operator * (const tensor& a, const tensor& b);
3655
friend tensor inv (const tensor& a, size_type d = 3);
3656
friend tensor diag (const point& d);
3658
// metric and geometric transformations:
3660
friend Float dotdot (const tensor&, const tensor&);
3661
friend Float norm2 (const tensor& a) { return dotdot(a,a); }
3662
friend Float dist2 (const tensor& a, const tensor& b) { return norm2(a-b); }
3663
friend Float norm (const tensor& a) { return ::sqrt(norm2(a)); }
3664
friend Float dist (const tensor& a, const tensor& b) { return norm(a-b); }
3665
Float determinant (size_type d = 3) const;
3666
friend Float determinant (const tensor& A, size_t d = 3);
3667
friend bool invert_3x3 (const tensor& A, tensor& result);
3670
// eigenvalues & eigenvectors:
3672
// a may be symmetric
3673
// where q=(q1,q2,q3) are eigenvectors in rows (othonormal matrix)
3674
// and d=(d1,d2,d3) are eigenvalues, sorted in decreasing order d1 >= d2 >= d3
3676
point eig (tensor& q, size_t dim = 3) const;
3677
point eig (size_t dim = 3) const;
3679
// singular value decomposition:
3681
// a can be unsymmetric
3682
// where u=(u1,u2,u3) are left pseudo-eigenvectors in rows (othonormal matrix)
3683
// v=(v1,v2,v3) are right pseudo-eigenvectors in rows (othonormal matrix)
3684
// and s=(s1,s2,s3) are eigenvalues, sorted in decreasing order s1 >= s2 >= s3
3686
point svd (tensor& u, tensor& v, size_t dim = 3) const;
3691
std::istream& get (std::istream&);
3694
tensor otimes (const point& a, const point& b, size_t na = 3);
3696
void cumul_otimes (tensor& t, const point& a, const point& b, size_t na = 3);
3697
void cumul_otimes (tensor& t, const point& a, const point& b, size_t na, size_t nb);
3700
File: rheolef.info, Node: point class, Up: Classes
3702
5.13 `point' - vertex of a mesh
3703
===============================
3705
(Source file: `nfem/basis/basic_point.h')
3710
Defines geometrical vertex as an array of coordinates. This array is
3711
also used as a vector of the three dimensional physical space.
3722
typedef size_t size_type;
3723
typedef T float_type;
3727
explicit basic_point () { _x[0] = T(); _x[1] = T(); _x[2] = T(); }
3729
explicit basic_point (
3733
{ _x[0] = x0; _x[1] = x1; _x[2] = x2; }
3736
basic_point<T>(const basic_point<T1>& p)
3737
{ _x[0] = p._x[0]; _x[1] = p._x[1]; _x[2] = p._x[2]; }
3740
basic_point<T>& operator = (const basic_point<T1>& p)
3741
{ _x[0] = p._x[0]; _x[1] = p._x[1]; _x[2] = p._x[2]; return *this; }
3745
T& operator[](int i_coord) { return _x[i_coord%3]; }
3746
const T& operator[](int i_coord) const { return _x[i_coord%3]; }
3747
T& operator()(int i_coord) { return _x[i_coord%3]; }
3748
const T& operator()(int i_coord) const { return _x[i_coord%3]; }
3753
std::istream& get (std::istream& s, int d = 3)
3756
case 1 : _x[1] = _x[2] = 0; return s >> _x[0];
3757
case 2 : _x[2] = 0; return s >> _x[0] >> _x[1];
3758
default: return s >> _x[0] >> _x[1] >> _x[2];
3762
std::ostream& put (std::ostream& s, int d = 3) const;
3764
// ccomparators: lexicographic order
3767
friend bool lexicographically_less (
3768
const basic_point<T>& a, const basic_point<T>& b) {
3769
for (size_type i = 0; i < d; i++) {
3770
if (a[i] < b[i]) return true;
3771
if (a[i] > b[i]) return false;
3773
return false; // equality
3777
friend bool operator == (const basic_point<T>& u, const basic_point<T>& v)
3778
{ return u[0] == v[0] && u[1] == v[1] && u[2] == v[2]; }
3780
basic_point<T>& operator+= (const basic_point<T>& v)
3781
{ _x[0] += v[0]; _x[1] += v[1]; _x[2] += v[2]; return *this; }
3783
basic_point<T>& operator-= (const basic_point<T>& v)
3784
{ _x[0] -= v[0]; _x[1] -= v[1]; _x[2] -= v[2]; return *this; }
3786
basic_point<T>& operator*= (const T& a)
3787
{ _x[0] *= a; _x[1] *= a; _x[2] *= a; return *this; }
3789
basic_point<T>& operator/= (const T& a)
3790
{ _x[0] /= a; _x[1] /= a; _x[2] /= a; return *this; }
3792
friend basic_point<T> operator+ (const basic_point<T>& u, const basic_point<T>& v)
3793
{ return basic_point<T> (u[0]+v[0], u[1]+v[1], u[2]+v[2]); }
3795
friend basic_point<T> operator- (const basic_point<T>& u)
3796
{ return basic_point<T> (-u[0], -u[1], -u[2]); }
3798
friend basic_point<T> operator- (const basic_point<T>& u, const basic_point<T>& v)
3799
{ return basic_point<T> (u[0]-v[0], u[1]-v[1], u[2]-v[2]); }
3802
friend basic_point<T> operator* (const T1& a, const basic_point<T>& u)
3803
{ return basic_point<T> (a*u[0], a*u[1], a*u[2]); }
3805
friend basic_point<T> operator* (const basic_point<T>& u, T a)
3806
{ return basic_point<T> (a*u[0], a*u[1], a*u[2]); }
3808
friend basic_point<T> operator/ (const basic_point<T>& u, const T& a)
3809
{ return basic_point<T> (u[0]/a, u[1]/a, u[2]/a); }
3811
friend basic_point<T> operator/ (const basic_point<T>& u, basic_point<T> v)
3812
{ return basic_point<T> (u[0]/v[0], u[1]/v[1], u[2]/v[2]); }
3814
friend basic_point<T> vect (const basic_point<T>& v, const basic_point<T>& w)
3815
{ return basic_point<T> (
3816
v[1]*w[2]-v[2]*w[1],
3817
v[2]*w[0]-v[0]*w[2],
3818
v[0]*w[1]-v[1]*w[0]); }
3821
// TODO: non-constant metric
3822
friend T dot (const basic_point<T>& u, const basic_point<T>& v)
3823
{ return u[0]*v[0]+u[1]*v[1]+u[2]*v[2]; }
3825
friend T norm2 (const basic_point<T>& u)
3826
{ return dot(u,u); }
3828
friend T norm (const basic_point<T>& u)
3829
{ return ::sqrt(norm2(u)); }
3831
friend T dist2 (const basic_point<T>& x, const basic_point<T>& y)
3832
{ return norm2(x-y); }
3834
friend T dist (const basic_point<T>& x, const basic_point<T>& y)
3835
{ return norm(x-y); }
3837
friend T dist_infty (const basic_point<T>& x, const basic_point<T>& y)
3838
{ return max(_my_abs(x[0]-y[0]),
3839
max(_my_abs(x[1]-y[1]), _my_abs(x[2]-y[2]))); }
3845
static T _my_abs(const T& x) { return (x > T(0)) ? x : -x; }
3848
T vect2d (const basic_point<T>& v, const basic_point<T>& w);
3851
T mixt (const basic_point<T>& u, const basic_point<T>& v, const basic_point<T>& w);
3853
// robust(exact) floating point predicates: return value as (0, > 0, < 0)
3854
// formally: orient2d(a,b,x) = vect2d(a-x,b-x)
3856
T orient2d(const basic_point<T>& a, const basic_point<T>& b,
3857
const basic_point<T>& x = basic_point<T>());
3859
// formally: orient3d(a,b,c,x) = mixt3d(a-x,b-x,c-x)
3861
T orient3d(const basic_point<T>& a, const basic_point<T>& b,
3862
const basic_point<T>& c, const basic_point<T>& x = basic_point<T>());
3865
File: rheolef.info, Node: ic0 class, Up: Classes
3867
5.14 `ic0' - incomplete Choleski factorization
3868
==============================================
3870
(Source file: `skit/lib/ic0.h')
3875
The class implements the icomplete Choleski factorization IC0(A) when
3876
$A$ is a square symmetric matrix stored in CSR format.
3877
csr<double> a = ...;
3878
ic0<double> fact(a);
3884
class basic_ic0 : public csr<T> {
3886
typedef typename csr<T>::size_type size_type;
3887
typedef T element_type;
3890
explicit basic_ic0 (const csr<T>& a);
3892
void inplace_solve (vec<T>& x) const;
3893
void solve (const vec<T>& b, vec<T>& x) const;
3894
vec<T> solve (const vec<T>& b) const;
3896
size_type nrow () const { return csr<T>::nrow(); }
3897
size_type ncol () const { return csr<T>::ncol(); }
3898
size_type nnz () const { return csr<T>::nnz(); }
3901
basic_ic0<T> ic0 (const csr<T>& a);
3904
File: rheolef.info, Node: permutation class, Up: Classes
3906
5.15 `permutation' - permutation matrix
3907
=======================================
3909
(Source file: `skit/lib/permutation.h')
3914
This class implements a permutation matrix. Permutation matrix are
3915
used in factorization implementation for optimizing the skyline format,
3916
as used by the implementation of the `ssk' class (*note ssk class::).
3918
Let a be a square invertible matrix in `csr' format (*note csr
3922
We get the optimal Gibbs permutation matrix using Gibbs, Pooles and
3923
Stockmeyer renumbering algorithm by:
3924
permutation p = gibbs(a);
3929
class permutation : public Array<std::vector<int>::size_type> {
3934
typedef Array<int>::size_type size_type;
3935
typedef Array<int>::difference_type difference_type;
3937
// allocators/deallocators:
3939
explicit permutation (size_type n = 0);
3942
permutation gibbs (const csr<T>& a);
3945
File: rheolef.info, Node: vec class, Up: Classes
3947
5.16 `vec' - dense vector
3948
=========================
3950
(Source file: `skit/lib/vec.h')
3955
The class implements an array. A declaration whithout any parametrers
3956
correspond to an array with a null size:
3958
Here, a `Float' array is specified. Note that the `Float' depends
3959
upon the configuration (*note Installing::). The constructor can be
3960
invocated whith a size parameter:
3962
Notes that the constructor can be invocated with an initializer:
3966
vec<Float> y = lambda;
3971
Linear combinations are `x+y', `x-y',
3972
`x*lambda', `lambda*x', `x/lambda'.
3974
Others combinations are `x*y', `x/y',
3977
The euclidian norm is `norm(x)' while `dot(x,y)'
3978
denotes the euclidian scalar product,
3980
Input/output routines overload "<<" and ">>" operators:
3984
*Note iorheo class::.
3989
Since the `vec<T>' class derives from the `Array<T>' class, the `vec'
3990
class present also a STL-like container interface.
3992
Actually, the `vec<T>' implementation uses a STL `vector<T>' class.
3998
class vec : public Array<T>
4002
typedef T element_type;
4003
typedef typename Array<T>::size_type size_type;
4005
// cstor, assignment
4006
explicit vec (unsigned int n = 0, const T& init_value = std::numeric_limits<T>::max())
4007
: Array<T>(n, init_value) {}
4008
vec<T> operator = (T lambda);
4011
unsigned int size () const { return Array<T>::size(); }
4012
unsigned int n () const { return size(); }
4014
#ifdef _RHEOLEF_HAVE_EXPRESSION_TEMPLATE
4017
vec(const VExpr<X>& x) : Array<T>(x.size()) { assign (*this, x); }
4020
vec<T>& operator = (const VExpr<X>& x)
4021
{ logmodif(*this); return assign (*this, x); }
4022
#endif // _RHEOLEF_HAVE_EXPRESSION_TEMPLATE
4025
template <class T> std::istream& operator >> (std::istream&, vec<T>&);
4026
template <class T> std::ostream& operator << (std::ostream&, const vec<T>&);
4029
File: rheolef.info, Node: ssk class, Up: Classes
4031
5.17 `ssk' - symmetric skyline matrix format
4032
============================================
4034
(Source file: `skit/lib/ssk.h')
4039
The class implements a symmetric matrix Choleski factorization. Let A
4040
be a square invertible matrix in `csr' format (*note csr class::).
4042
We get the factorization by:
4043
ssk<Float> m = ldlt(a);
4044
Each call to the direct solver for a*x = b writes either:
4045
vec<Float> x = m.solve(b);
4052
The storage is either skyline, multi-file or chevron. This
4053
alternative depends upon the configuration (*note Installing::). The
4054
chevron data structure is related to the multifrontal algorithm and is
4055
implemented by using the spooles library. The multi-file data
4056
structure refers to the out-of-core algorithm and is implemented by
4057
using the taucs library. If such a library is not available, the
4058
`ssk' class implements a skyline storage scheme and the profil storage
4059
of the `ssk' matrix is optimized by optimizing the renumbering, by
4060
using the Gibbs, Pooles and Stockmeyer algorithm available via the
4061
`permutation' class (*note permutation class::).
4063
Algorithm 582, collected Algorithms from ACM. Algorithm
4064
appeared in ACM-Trans. Math. Software, vol.8, no. 2, Jun.,
4066
When implementing the skyline data structure, we can go back to the
4067
`csr' format, for the pupose of pretty-printing:
4068
cout << ps << color << logscale << csr<Float>(m);
4074
class ssk : smart_pointer<spooles_rep<T> > {
4078
typedef T element_type;
4079
typedef typename spooles_rep<T>::size_type size_type;
4081
// allocators/deallocators:
4084
explicit ssk<T> (const csr<T>&);
4088
size_type nrow () const;
4089
size_type ncol () const;
4090
size_type nnz () const;
4092
// factorisation and solver:
4094
void solve(const vec<T>& b, vec<T>& x) const;
4095
vec<T> solve(const vec<T>& b) const;
4096
void factorize_ldlt();
4098
const spooles_rep<T>& data() const {
4099
return smart_pointer<spooles_rep<T> >::data();
4101
spooles_rep<T>& data() {
4102
return smart_pointer<spooles_rep<T> >::data();
4106
ssk<T> ldlt (const csr<T>& m);
4112
class ssk : smart_pointer<umfpack_rep<T> > {
4116
typedef T element_type;
4117
typedef typename umfpack_rep<T>::size_type size_type;
4119
// allocators/deallocators:
4122
explicit ssk<T> (const csr<T>&);
4126
size_type nrow () const;
4127
size_type ncol () const;
4128
size_type nnz () const;
4130
// factorisation and solver:
4132
void solve(const vec<T>& b, vec<T>& x) const;
4133
vec<T> solve(const vec<T>& b) const;
4134
void factorize_ldlt();
4135
void factorize_lu();
4137
const umfpack_rep<T>& data() const {
4138
return smart_pointer<umfpack_rep<T> >::data();
4140
umfpack_rep<T>& data() {
4141
return smart_pointer<umfpack_rep<T> >::data();
4145
ssk<T> ldlt (const csr<T>& m);
4148
File: rheolef.info, Node: csr class, Up: Classes
4150
5.18 `csr' - compressed sparse row matrix format
4151
================================================
4153
(Source file: `skit/lib/csr.h')
4158
The class implements a matrix in compressed sparse row format. A
4159
declaration whithout any parametrers correspond to a matrix with null
4162
Notes that the constructor can be invocated with an initializer:
4165
Input and output, as usual (*note iorheo class::):
4169
Default is the Harwell-Boeing format. Various others formated options
4170
are available: matlab and postscript plots.
4172
Affectation from a scalar
4175
The matrix/vector multiply:
4177
and the transposed matrix/ vector multiply:
4180
The binary operators are:
4181
a*b, a+b, a-b, lambda*a, a*lambda, a/lambda
4182
The unary operators are sign inversion and transposition:
4184
The combination with a diagonal matrix is not yet completely
4185
available. The interface would be something like:
4186
basic_diag<double> d;
4190
left_div(a,d) // aij/dii
4191
When applied to the matrix directly: this feature is not yet
4192
completely available. The interface would be something like:
4196
a.left_mult(d); // a := d*a
4197
a /= d; // aij := aij/djj
4198
a.left_div(d); // aij := aij/dii
4199
The combination with a constant-identity matrix: this feature is not
4200
yet completely available. The interface would be something like:
4202
a + k*EYE, k*EYE + a, a - k*EYE, k*EYE - a,
4206
Get the lower triangular part:
4207
csr<double> l = tril(a);
4208
Conversely, `triu' get the upper triangular part.
4210
For optimizing the profile storage, I could be convenient to use a
4211
renumbering algorithm (see also *note ssk class:: and *note
4212
permutation class::).
4213
permutation p = gibbs(a);
4214
csr<double> b = perm(a, p, q); // B := P*A*trans(Q)
4215
Horizontal matrix concatenation:
4217
Vertical matrix concatenation:
4219
Explicit conversion from an associative `asr e' matrix writes:
4221
from a dense `dns m' matrix writes:
4228
The `csr' class is currently under reconstruction for the distributed
4229
memory support by using a MPI-based implementation.
4232
File: rheolef.info, Node: diag class, Up: Classes
4234
5.19 `basic_diag' - diagonal matrix
4235
===================================
4237
(Source file: `skit/lib/diag.h')
4242
The class implements a diagonal matrix. A declaration whithout any
4243
parametrers correspond to a null size matrix:
4244
basic_diag<Float> d;
4245
The constructor can be invocated whith a size parameter:
4246
basic_diag<Float> d(n);
4247
or an initialiser, either a vector (*note vec class::):
4248
basic_diag<Float> d = basic_diag(v);
4249
or a csr matrix *note csr class:::
4250
basic_diag<Float> d = basic_diag(a);
4251
The conversion from `basic_diag' to `vec' or `csr' is explicit.
4253
When a diagonal matrix is constructed from a `csr' matrix, the
4254
definition of the diagonal of matrix is _always_ a vector of size NROW
4255
which contains the elements in rows 1 to NROW of the matrix that are
4256
contained in the diagonal. If the diagonal element falls outside the
4257
matrix, i.e. NCOL < NROW then it is defined as a zero entry.
4262
Since the `basic_diag' class derives from the `vec', the `basic_diag'
4263
class present also a STL-like interface.
4269
class basic_diag : public vec<T> {
4274
typedef typename vec<T>::element_type element_type;
4275
typedef typename vec<T>::size_type size_type;
4276
typedef typename vec<T>::iterator iterator;
4278
// allocators/deallocators:
4280
explicit basic_diag (size_type sz = 0);
4281
explicit basic_diag (const vec<T>& u);
4282
explicit basic_diag (const csr<T>& a);
4286
basic_diag<T> operator = (const T& lambda);
4290
size_type nrow () const { return vec<T>::size(); }
4291
size_type ncol () const { return vec<T>::size(); }
4293
// basic_diag as a preconditionner: solves D.x=b
4295
vec<T> solve (const vec<T>& b) const;
4296
vec<T> trans_solve (const vec<T>& b) const;
4299
basic_diag<T> dcat (const basic_diag<T>& a1, const basic_diag<T>& a2);
4302
basic_diag<T> operator / (const T& lambda, const basic_diag<T>& d);
4306
operator * (const basic_diag<T>& d, const vec<T>& x);
4309
vec<T> left_div (const vec<T>& x, const basic_diag<T>& d);
4312
File: rheolef.info, Node: iorheo class, Up: Classes
4314
5.20 `iorheo' - input and output functions and manipulation
4315
===========================================================
4317
(Source file: `util/lib/iorheo.h')
4322
input geo in standard file format:
4324
output geo in standard file format:
4326
output geo in gnuplot format:
4327
cout << gnuplot << g;
4332
output manipulators enable the selection of some pretty graphic
4333
options, in an elegant fashion.
4335
Boolean manipulators
4336
--------------------
4338
The boolean manipulators set an internal optional flag. Each
4339
option has its negative counterpart, as `verbose' and `noverbose',
4340
by adding the `no' prefix.
4341
cout << noverbose << a;
4344
trace some details, such as loading, storing or unix
4345
command operations on `cerr'. Default is on.
4348
delete temporary files during graphic outputs. Default is
4352
run unix operations, such as `gnuplot' or `plotmtv' or
4353
`vtk'. Note that the corresponding files are created.
4357
perform transposition when loading/soring a `csr' matrix
4358
from Harwell-Boeing file. This feature is available, since
4359
the file format store matrix in transposed format.
4363
when using matrix sparse postscript plot manipulator `ps'
4364
and `color'. The color scale is related to a logarithmic
4377
when using the `vtk' or `mayavi' manipulators for a mesh or
4381
volume rendering by using ray cast functions.
4385
Vector-valued fields are rendered by using arrows (`velocity')
4386
or deformed meshes (`deformation'). For `vtk' or
4387
`plotmtv' rendering.
4390
Scalar valued fields in two dimension are rendered by using
4391
a tridimensionnal surface elevation. For `vtk' or
4392
`plotmtv' rendering.
4396
try to reuse the supplied space. Default is on.
4398
File format manipulators
4399
------------------------
4401
The "format" manipulator group applies for streams. Its value is
4402
an enumerated type, containing the following possibilities:
4404
use the default textual input/output file format. For
4405
instance, this is `.geo' for meshes, `.field' for
4406
discretized functions. This default format is specified
4407
in the corresponding class documentation (see also
4408
*note geo class:: and *note field class::). This is the
4412
uses `.bamg' Frederic Hecht's bidimensional anisotropic
4413
mesh generator file format for geo input/output operation.
4416
uses `.node' `.ele' and `.face' Hang Si's tridimensional
4417
mesh generator file format for geo input/output operation.
4420
uses `.mmg3d' Cecile Dobrzynski's tridimensional anisotropic
4421
mesh generator file format for geo input/output operation.
4424
uses `.msh' gmsh Christophe Geuzaine and Jean-François
4425
Remacle mesh generator file format for geo input/output
4429
uses `.m' (bidimensional) or `.v' (tridimensionnal) Carl
4430
Ollivier-Gooch 's mesh generator file format for geo input/output
4434
uses `.qmg' Stephen A. Vavasis's mesh generator file
4435
format for geo input/output operation.
4438
uses `.vtk' mesh file format for geo input/output
4439
operations. This file format is suitable for graphic
4443
uses `.vtk' polydata (specific for polygonal boundaries)
4444
mesh file format for geo input/output operations. This
4445
file format is suitable for graphic treatment.
4448
uses `.cemagref' surface mesh (topography, with a `z' cote).
4449
This file format is used at Cemagref (french research center for
4450
mountains, `http://www.cemagref.fr').
4453
output an extensive listing of the class data structure.
4454
This option is used for debugging purpose.
4457
uses `.hb' Harwell-Boeing file format for sparse matrix
4458
input/output operation. This is the default.
4461
uses `.mm' Matrix-Market file format for sparse matrix input/output
4466
uses `.m' Matlab file format for sparse matrix output
4470
uses `.m' Matlab sparse file format for sparse matrix output
4474
uses `.ps' postscript for sparse matrix output operation.
4477
for mesh and field outputs. Generate `.vtk' data file and
4478
the `.tcl' command script file and run the `vtk'
4479
command on the `.tcl'.
4482
for field outputs. Generate `.vtk' data file and the
4483
`.py' command script file and run the `python' command
4484
on the `.py' associated to the `mayavi'/`vtk' library.
4487
for boundary cad outputs. Generate `.off' data file
4488
and run the `geomview' command.
4491
for mesh and field outputs. Generate `.gdat' data file
4492
and the `.plot' command script file and run the
4493
`gnuplot' command on the `.plot'.
4496
for mesh and field outputs. Generate `.mtv' data file
4497
and run the `plotmtv' command.
4500
for mesh output. Generate `.x3d' data file and
4501
run the `x3d' command. This tool has fast rotation
4505
for mesh output. Generate `.atom' data file and
4506
run the `PlotM' command. Tridimensional mesh rendering is
4507
performed as a chemical molecule: nodes as balls and
4508
edges as tubes. The `PlotM' tool is developped
4509
at Cornell University Laboratory of Atomic and Solid State
4510
Physics (LASSP) in a Joint Study with IBM, with support by
4511
the Materials Science Center and Corning Glassworks.
4516
The "color" manipulator group acts for sparse matrix postscript
4517
output. Its value is an enumerated type, containing three
4521
cout << black_and_white << c;
4522
The default is to generate a color postcript file. Conversely,
4523
its act for field rendering, via `mayavi'.
4525
Valuated manipulators
4526
---------------------
4528
Some manipulators takes an agument that specifies a value.
4529
cout << geomview << bezieradapt << subdivide(5) << my_cad_boundary;
4530
cout << vtk << iso << isovalue(2.5) << my_mesh;
4531
cout << velocity << plotmtv << vectorscale(0.1) << uh;
4532
See also *note catchmark class:: for input-output of vector-valued
4536
`n_isovalue_negative INT'
4541
File: rheolef.info, Node: Vector class, Up: Classes
4543
5.21 `Vector' - STL `vector<T>' with reference counting
4544
=======================================================
4546
(Source file: `util/lib/Vector.h')
4551
The class implement a reference counting wrapper for the STL
4552
`vector<T>' container class, with shallow copies. See also:
4554
The Standard Template Library, by Alexander Stephanov and Meng
4557
This class provides the full `vector<T>' interface specification an
4558
could be used instead of `vector<T>'.
4564
T& operator[](size_type)
4565
as in `v[i]' may checks the reference count for each access. For a
4566
loop, a better usage is:
4567
Vector<T>::iterator i = v.begin();
4568
Vector<T>::iterator last = v.end();
4569
while (i != last) { ...}
4570
and the reference count check step occurs only two time, when
4571
accessing via `begin()' and `end()'.
4573
Thus, in order to encourage users to do it, we declare private theses
4574
member functions. A synonym of `operator[]' is `at'.
4580
class Vector : private smart_pointer<vector_rep<T> > {
4587
typedef const_iterator;
4590
typedef const_reference;
4592
typedef difference_type;
4594
typedef reverse_iterator;
4595
typedef const_reverse_iterator;
4597
// allocation/deallocation:
4599
explicit Vector (size_type n = 0, const T& value = T ());
4600
Vector (const_iterator first, const_iterator last);
4601
void reserve (size_type n);
4602
void swap (Vector<T>& x) ;
4607
const_iterator begin () const;
4609
const_iterator end () const;
4610
reverse_iterator rbegin();
4611
const_reverse_iterator rbegin() const;
4612
reverse_iterator rend();
4613
const_reverse_iterator rend() const;
4614
size_type size () const;
4615
size_type max_size () const;
4616
size_type capacity () const;
4617
bool empty () const;
4618
void resize (size_type sz, T v = T ()); // non-standard ?
4620
const_reference operator[] (size_type n) const;
4621
reference operator[] (size_type n);
4623
const_reference at (size_type n) const; // non-standard ?
4624
reference at (size_type n);
4626
const_reference front () const;
4628
const_reference back () const;
4632
void push_back (const T& x);
4633
iterator insert (iterator position, const T& x = T ());
4634
void insert (iterator position, size_type n, const T& x);
4635
void insert (iterator position, const_iterator first, const_iterator last);
4637
void erase (iterator position);
4638
void erase (iterator first, iterator last);
4642
File: rheolef.info, Node: catchmark class, Up: Classes
4644
5.22 `catchmark' - iostream manipulator
4645
=======================================
4647
(Source file: `util/lib/catchmark.h')
4652
The `catchmark' is used for building labels used for input-output
4653
of vector-valued fields (see *note field class::):
4654
cin >> catchmark("f") >> fh;
4655
cout << catchmark("u") << uh
4656
<< catchmark("w") << wh
4657
<< catchmark("psi") << psih;
4658
Assuming its value for output is `"u"', the corresponding
4659
labels will be `"#u0"', `"#u1"', `"#u2"', ...
4666
catchmark(const std::string& x);
4667
friend std::istream& operator >> (std::istream& is, const catchmark& m);
4668
friend std::ostream& operator << (std::ostream& os, const catchmark& m);
4674
File: rheolef.info, Node: rheostream class, Up: Classes
4676
5.23 `irheostream', `orheostream' - large data streams
4677
======================================================
4679
(Source file: `util/lib/rheostream.h')
4684
This class provides a stream interface for large data management.
4685
File decompresion is assumed using `gzip' and a recursive seach in a
4686
directory list is provided for input.
4687
orheostream foo("NAME", "suffix");
4689
ofstream foo("NAME.suffix").
4690
However, if NAME does not end with `.suffix', then `.suffix' is
4691
added, and compression is done with gzip, adding an additional `.gz'
4695
irheostream foo("NAME","suffix");
4697
ifstream foo("NAME.suffix").
4698
However, we look at a search path environment variable `RHEOPATH' in
4699
order to find NAME while suffix is assumed. Moreover, `gzip'
4700
compressed files, ending with the `.gz' suffix is assumed, and
4701
decompression is done.
4703
Finally, a set of useful functions are provided.
4709
irheostream is("results", "data");
4710
will recursively look for a `results[.data[.gz]]' file in the
4711
directory mentionned by the `RHEOPATH' environment variable.
4713
For instance, if you insert in our ".cshrc" something like:
4714
setenv RHEOPATH ".:/home/dupont:/usr/local/math/demo"
4715
the process will study the current directory `.', then, if neither
4716
`square.data.gz' nor `square.data' exits, it scan all subdirectory of
4717
the current directory. Then, if file is not founded, it start
4718
recusively in `/home/dupond' and then in `/usr/local/math/demo'.
4720
File decompression is performed by using the `gzip' command, and data
4721
are pipe-lined directly in memory.
4723
If the file start with `.' as `./square' or with a `/' as
4724
`/home/oscar/square', no search occurs and `RHEOPATH' environment
4725
variable is not used.
4727
Also, if the environment variable `RHEOPATH' is not set, the default
4728
value is the current directory `.'.
4731
orheostream os("newresults", "data");
4732
file compression is assumed, and "newresults.data.gz" will be created.
4734
File loading and storing are mentionned by a message, either:
4735
! load "./results.data.gz"
4737
! file "./newresults.data.gz" created.
4738
on the `clog' stream. By adding the following:
4740
you turn off these messages (*note iorheo class::).
4745
class irheostream : public io::filtering_stream<io::input> {
4747
irheostream() : io::filtering_stream<io::input>() {}
4748
irheostream(const std::string& name, const std::string& suffix = std::string());
4749
virtual ~irheostream();
4750
void open (const std::string& name, const std::string& suffix = std::string());
4755
class orheostream : public io::filtering_stream<io::output> {
4757
orheostream() : io::filtering_stream<io::output>() {}
4758
orheostream(const std::string& name, const std::string& suffix = std::string());
4759
virtual ~orheostream();
4760
void open (const std::string& name, const std::string& suffix = std::string());
4765
std::string itos (std::string::size_type i);
4766
std::string ftos (const Float& x);
4768
// catch first occurence of string in file
4769
bool scatch (std::istream& in, const std::string& ch);
4771
// has_suffix("toto.suffix", "suffix") -> true
4772
bool has_suffix (const std::string& name, const std::string& suffix);
4774
// "toto.suffix" --> "toto"
4775
std::string delete_suffix (const std::string& name, const std::string& suffix);
4777
// "/usr/local/dir/toto.suffix" --> "toto.suffix"
4778
std::string get_basename (const std::string& name);
4780
// "/usr/local/dir/toto.suffix" --> "/usr/local/dir"
4781
std::string get_dirname (const std::string& name);
4783
// "toto" --> "/usr/local/math/data/toto.suffix"
4784
std::string get_full_name_from_rheo_path (const std::string& rootname, const std::string& suffix);
4786
// "." + "../geodir" --> ".:../geodir"
4787
void append_dir_to_rheo_path (const std::string& dir);
4789
// "../geodir" + "." --> "../geodir:."
4790
void prepend_dir_to_rheo_path (const std::string& dir);
4792
bool file_exists (const std::string& filename);
4795
bool is_float (const std::string&);
4796
Float to_float (const std::string&);
4799
File: rheolef.info, Node: Algorithms, Up: Top
4806
* newton algorithm::
4807
* damped-newton algorithm::
4808
* riesz_representer algorithm::
4809
* mixed_solver algorithm::
4810
* pminres algorithm::
4813
* bicgstab algorithm::
4814
* puzawa algorithm::
4818
File: rheolef.info, Node: newton algorithm, Up: Algorithms
4820
6.1 `newton' - Newton nonlinear algorithm
4821
=========================================
4823
(Source file: `nfem/lib/newton.h')
4828
Nonlinear Newton algorithm for the resolution of the following problem:
4830
A simple call to the algorithm writes:
4833
newton (P, uh, tol, max_iter);
4834
The `my_problem' class may contains methods for the evaluation of F
4835
(aka residue) and its derivative:
4839
field residue (const field& uh) const;
4840
void update_derivative (const field& uh) const;
4841
field derivative_solve (const field& mrh) const;
4842
Float norm (const field& uh) const;
4843
Float dual_norm (const field& Muh) const;
4845
See the example `p-laplacian.h' in the user's documentation for
4851
template <class Problem, class Field>
4852
int newton (Problem P, Field& uh, Float& tol, size_t& max_iter, std::ostream *p_cerr = 0) {
4853
if (p_cerr) *p_cerr << "# Newton: n r" << std::endl;
4854
for (size_t n = 0; true; n++) {
4855
Field rh = P.residue(uh);
4856
Float r = P.dual_norm(rh);
4857
if (p_cerr) *p_cerr << n << " " << r << std::endl;
4858
if (r <= tol) { tol = r; max_iter = n; return 0; }
4859
if (n == max_iter) { tol = r; return 1; }
4860
P.update_derivative (uh);
4861
Field delta_uh = P.derivative_solve (-rh);
4867
File: rheolef.info, Node: damped-newton algorithm, Up: Algorithms
4869
6.2 `damped_newton' - damped Newton nonlinear algorithm
4870
=======================================================
4872
(Source file: `nfem/lib/damped-newton.h')
4877
Nonlinear damped Newton algorithm for the resolution of the following
4880
A simple call to the algorithm writes:
4883
damped_newton (P, uh, tol, max_iter);
4884
The `my_problem' class may contains methods for the evaluation of F
4885
(aka residue) and its derivative:
4889
field residue (const field& uh) const;
4890
void update_derivative (const field& uh) const;
4891
field derivative_trans_mult (const field& mrh) const;
4892
field derivative_solve (const field& mrh) const;
4893
Float norm (const field& uh) const;
4894
Float dual_norm (const field& Muh) const;
4896
See the example `p-laplacian.h' in the user's documentation for
4902
template <class Problem, class Field, class Real, class Size>
4903
int damped_newton (Problem P, Field& u, Real& tol, Size& max_iter, std::ostream* p_cerr=0) {
4904
return damped_newton(P, newton_identity_preconditioner(), u, tol, max_iter, p_cerr);
4908
File: rheolef.info, Node: riesz_representer algorithm, Up: Algorithms
4910
6.3 `riesz_representer' - integrate a function by using quadrature formulae
4911
===========================================================================
4913
(Source file: `nfem/lib/riesz_representer.h')
4918
The function `riesz_representer' implements the approximation of an
4919
integral by using quadrature formulae. This is an expreimental
4920
implementation: please do not use yet for practical usage.
4925
template <class Function> field riesz_representer (const space& Vh,
4928
template <class Function> field riesz_representer (const space& Vh,
4929
const Function& f, quadrature_option_type qopt);
4934
The following code compute the Riesz representant, denoted by `mfh' of
4935
f(x), and the integral of f over the domain omega:
4936
Float f(const point& x);
4938
space Vh (omega_h, "P1");
4939
field mfh = riesz_representer(Vh, f);
4940
Float int_f = dot(mfh, field(Vh,1.0));
4941
The Riesz representer is the `mfh' vector of values:
4942
mfh(i) = integrate f(x) phi_i(x) dx
4943
where phi_i is the i-th basis function in Vh and the integral is
4944
evaluated by using a quadrature formulae. By default the quadrature
4945
formule is the Gauss one with the order equal to the polynomial order
4946
of Vh. Alternative quadrature formulae and order is available by
4947
passing an optional variable to riesz_representer.
4952
template <class Function>
4957
quadrature_option_type qopt
4958
= quadrature_option_type(quadrature_option_type::max_family,0))
4961
File: rheolef.info, Node: mixed_solver algorithm, Up: Algorithms
4963
6.4 `pcg_abtb', `pcg_abtbc', `pminres_abtb', `pminres_abtbc' - solvers for mixed linear problems
4964
================================================================================================
4966
(Source file: `skit/lib/mixed_solver.h')
4971
template <class Matrix, class Vector, class Solver, class Preconditioner, class Size, class Real>
4972
int pcg_abtb (const Matrix& A, const Matrix& B, Vector& u, Vector& p,
4973
const Vector& Mf, const Vector& Mg, const Preconditioner& S1,
4974
const Solver& inner_solver_A, Size& max_iter, Real& tol,
4975
std::ostream *p_cerr = 0, std::string label = "pcg_abtb");
4977
template <class Matrix, class Vector, class Solver, class Preconditioner, class Size, class Real>
4978
int pcg_abtbc (const Matrix& A, const Matrix& B, const Matrix& C, Vector& u, Vector& p,
4979
const Vector& Mf, const Vector& Mg, const Preconditioner& S1,
4980
const Solver& inner_solver_A, Size& max_iter, Real& tol,
4981
std::ostream *p_cerr = 0, std::string label = "pcg_abtbc");
4982
The synopsis is the same with the pminres algorithm.
4987
See the user's manual for practical examples for the nearly
4988
incompressible elasticity, the Stokes and the Navier-Stokes problems.
4993
Preconditioned conjugate gradient algorithm on the pressure p applied to
4994
the stabilized stokes problem:
4995
[ A B^T ] [ u ] [ Mf ]
4997
[ B -C ] [ p ] [ Mg ]
4998
where A is symmetric positive definite and C is symmetric positive
4999
and semi-definite. Such mixed linear problems appears for instance
5000
with the discretization of Stokes problems with stabilized P1-P1
5001
element, or with nearly incompressible elasticity. Formaly u =
5002
inv(A)*(Mf - B^T*p) and the reduced system writes for all
5003
non-singular matrix S1:
5004
inv(S1)*(B*inv(A)*B^T)*p = inv(S1)*(B*inv(A)*Mf - Mg)
5005
Uzawa or conjugate gradient algorithms are considered on the
5006
reduced problem. Here, S1 is some preconditioner for the Schur
5007
complement S=B*inv(A)*B^T. Both direct or iterative solvers for S1*q
5008
= t are supported. Application of inv(A) is performed via a call to
5009
a solver for systems such as A*v = b. This last system may be
5010
solved either by direct or iterative algorithms, thus, a general
5011
matrix solver class is submitted to the algorithm. For most
5012
applications, such as the Stokes problem, the mass matrix for the p
5013
variable is a good S1 preconditioner for the Schur complement. The
5014
stoping criteria is expressed using the S1 matrix, i.e. in L2 norm
5015
when this choice is considered. It is scaled by the L2 norm of the
5016
right-hand side of the reduced system, also in S1 norm.
5019
File: rheolef.info, Node: pminres algorithm, Up: Algorithms
5021
6.5 `pminres' - conjugate gradient algorithm.
5022
=============================================
5024
(Source file: `skit/lib/pminres.h')
5029
template <class Matrix, class Vector, class Preconditioner, class Real>
5030
int pminres (const Matrix &A, Vector &x, const Vector &b, const Preconditioner &M,
5031
int &max_iter, Real &tol, std::ostream *p_cerr=0);
5036
The simplest call to 'pminres' has the folling form:
5037
size_t max_iter = 100;
5039
int status = pminres(a, x, b, EYE, max_iter, tol, &cerr);
5044
`pminres' solves the symmetric positive definite linear system Ax=b
5045
using the Conjugate Gradient method.
5047
The return value indicates convergence within max_iter (input)
5048
iterations (0), or no convergence within max_iter iterations (1).
5049
Upon successful return, output arguments have the following values:
5051
approximate solution to Ax = b
5054
the number of iterations performed before the tolerance was reached
5057
the residual after the final iteration
5062
`pminres' follows the algorithm described in "Solution of sparse
5063
indefinite systems of linear equations", C. C. Paige and M. A.
5064
Saunders, SIAM J. Numer. Anal., 12(4), 1975. For more, see
5065
http://www.stanford.edu/group/SOL/software.html and also the PhD
5066
"Iterative methods for singular linear equations and least-squares
5067
problems", S.-C. T. Choi, Stanford University, 2006,
5068
http://www.stanford.edu/group/SOL/dissertations/sou-cheng-choi-thesis.pdf
5069
at page 60. The present implementation style is inspired from `IML++
5070
1.2' iterative method library, `http://math.nist.gov/iml++'.
5075
template <class Matrix, class Vector, class Preconditioner, class Real, class Size>
5076
int pminres(const Matrix &A, Vector &x, const Vector &Mb, const Preconditioner &M,
5077
Size &max_iter, Real &tol, std::ostream *p_cerr = 0, std::string label = "minres")
5079
Vector b = M.solve(Mb);
5080
Real norm_b = sqrt(fabs(dot(Mb,b)));
5081
if (norm_b == Real(0.)) norm_b = 1;
5083
Vector Mr = Mb - A*x;
5084
Vector z = M.solve(Mr);
5085
Real beta2 = dot(Mr, z);
5086
Real norm_r = sqrt(fabs(beta2));
5087
if (p_cerr) (*p_cerr) << "[" << label << "] #iteration residue" << std::endl;
5088
if (p_cerr) (*p_cerr) << "[" << label << "] 0 " << norm_r/norm_b << std::endl;
5089
if (beta2 < 0 || norm_r <= tol*norm_b) {
5090
tol = norm_r/norm_b;
5092
warning_macro ("beta2 = " << beta2 << " < 0: stop");
5095
Real beta = sqrt(beta2);
5097
Vector Mv = Mr/beta;
5103
Vector u_old (x.size(), 0.);
5104
Vector Mv_old (x.size(), 0.);
5105
Vector w (x.size(), 0.);
5106
Vector w_old (x.size(), 0.);
5107
Vector w_old2 (x.size(), 0.);
5108
for (Size n = 1; n <= max_iter; n++) {
5112
Real alpha = dot(Mr, u);
5113
Mr = Mr - alpha*Mv - beta*Mv_old;
5114
z = z - alpha*u - beta*u_old;
5117
warning_macro ("pminres: machine precision problem");
5118
tol = norm_r/norm_b;
5122
Real beta_old = beta;
5125
Real c_old2 = c_old;
5126
Real s_old2 = s_old;
5129
Real rho0 = c_old*alpha - c_old2*s_old*beta_old;
5130
Real rho2 = s_old*alpha + c_old2*c_old*beta_old;
5131
Real rho1 = sqrt(sqr(rho0) + sqr(beta));
5132
Real rho3 = s_old2 * beta_old;
5139
w = (u - rho2*w_old - rho3*w_old2)/rho1;
5148
if (p_cerr) (*p_cerr) << "[" << label << "] " << n << " " << norm_r/norm_b << std::endl;
5149
if (norm_r <= tol*norm_b) {
5150
tol = norm_r/norm_b;
5155
tol = norm_r/norm_b;
5160
File: rheolef.info, Node: qmr algorithm, Up: Algorithms
5162
6.6 `qmr' - quasi-minimal residual algoritm
5163
===========================================
5165
(Source file: `skit/lib/qmr.h')
5170
template <class Matrix, class Vector, class Preconditioner1,
5171
class Preconditioner2, class Real>
5172
int qmr (const Matrix &A, Vector &x, const Vector &b,
5173
const Preconditioner1 &M1, const Preconditioner2 &M2,
5174
int &max_iter, Real &tol);
5179
The simplest call to 'qmr' has the folling form:
5180
int status = qmr(a, x, b, EYE, EYE, 100, 1e-7);
5185
`qmr' solves the unsymmetric linear system Ax = b using the the
5186
quasi-minimal residual method.
5188
The return value indicates convergence within max_iter (input)
5189
iterations (0), or no convergence within max_iter iterations (1).
5190
Upon successful return, output arguments have the following values:
5192
approximate solution to Ax = b
5195
the number of iterations performed before the tolerance was reached
5198
the residual after the final iteration
5200
A return value of 1 indicates that the method did not reach the
5201
specified convergence tolerance in the maximum numbefr of iterations.
5202
A return value of 2 indicates that a breackdown associated with `rho'
5203
occurred. A return value of 3 indicates that a breackdown associated
5204
with `beta' occurred. A return value of 4 indicates that a
5205
breackdown associated with `gamma' occurred. A return value of 5
5206
indicates that a breackdown associated with `delta' occurred. A
5207
return value of 6 indicates that a breackdown associated with `epsilon'
5208
occurred. A return value of 7 indicates that a breackdown associated
5214
`qmr' is an iterative template routine.
5216
`qmr' follows the algorithm described on p. 24 in
5218
Templates for the Solution of Linear Systems: Building
5219
Blocks for Iterative Methods, 2nd Edition, R.
5220
Barrett, M. Berry, T. F. Chan, J. Demmel, J. Donato, J. Dongarra,
5221
V. Eijkhout, R. Pozo, C. Romine, H. Van der Vorst,
5222
SIAM, 1994, `ftp.netlib.org/templates/templates.ps'.
5224
The present implementation is inspired from `IML++ 1.2' iterative
5225
method library, `http://math.nist.gov/iml++'.
5228
File: rheolef.info, Node: pcg algorithm, Up: Algorithms
5230
6.7 `pcg' - conjugate gradient algorithm.
5231
=========================================
5233
(Source file: `skit/lib/pcg.h')
5238
template <class Matrix, class Vector, class Preconditioner, class Real>
5239
int pcg (const Matrix &A, Vector &x, const Vector &b,
5240
const Preconditioner &M, int &max_iter, Real &tol, std::ostream *p_cerr=0);
5245
The simplest call to 'pcg' has the folling form:
5246
size_t max_iter = 100;
5248
int status = pcg(a, x, b, EYE, max_iter, tol, &cerr);
5253
`pcg' solves the symmetric positive definite linear system Ax=b using
5254
the Conjugate Gradient method.
5256
The return value indicates convergence within max_iter (input)
5257
iterations (0), or no convergence within max_iter iterations (1).
5258
Upon successful return, output arguments have the following values:
5260
approximate solution to Ax = b
5263
the number of iterations performed before the tolerance was reached
5266
the residual after the final iteration
5271
`pcg' is an iterative template routine.
5273
`pcg' follows the algorithm described on p. 15 in
5275
Templates for the Solution of Linear Systems: Building
5276
Blocks for Iterative Methods, 2nd Edition, R.
5277
Barrett, M. Berry, T. F. Chan, J. Demmel, J. Donato, J. Dongarra,
5278
V. Eijkhout, R. Pozo, C. Romine, H. Van der Vorst,
5279
SIAM, 1994, `ftp.netlib.org/templates/templates.ps'.
5281
The present implementation is inspired from `IML++ 1.2' iterative
5282
method library, `http://math.nist.gov/iml++'.
5287
template < class Matrix, class Vector, class Preconditioner, class Real, class Size>
5288
int pcg(const Matrix &A, Vector &x, const Vector &Mb, const Preconditioner &M,
5289
Size &max_iter, Real &tol, std::ostream *p_cerr = 0, std::string label = "cg")
5291
Vector b = M.solve(Mb);
5292
Real norm2_b = dot(Mb,b);
5293
if (norm2_b == Real(0)) norm2_b = 1;
5294
Vector Mr = Mb - A*x;
5296
if (p_cerr) (*p_cerr) << "[" << label << "] #iteration residue" << std::endl;
5298
for (Size n = 0; n <= max_iter; n++) {
5299
Vector r = M.solve(Mr);
5300
Real prev_norm2_r = norm2_r;
5301
norm2_r = dot(Mr, r);
5302
if (p_cerr) (*p_cerr) << "[" << label << "] " << n << " " << ::sqrt(norm2_r/norm2_b) << std::endl;
5303
if (norm2_r <= sqr(tol)*norm2_b) {
5304
tol = ::sqrt(norm2_r/norm2_b);
5311
Real beta = norm2_r/prev_norm2_r;
5315
Real alpha = norm2_r/dot(Mq, p);
5319
tol = ::sqrt(norm2_r/norm2_b);
5324
File: rheolef.info, Node: bicgstab algorithm, Up: Algorithms
5326
6.8 `bicgstab' - bi-conjugate gradient stabilized method
5327
========================================================
5329
(Source file: `skit/lib/bicgstab.h')
5334
int bicgstab (const Matrix &A, Vector &x, const Vector &b,
5335
const Preconditioner &M, int &max_iter, Real &tol);
5340
The simplest call to 'bicgstab' has the folling form:
5341
int status = bicgstab(a, x, b, EYE, 100, 1e-7);
5346
`bicgstab' solves the unsymmetric linear system Ax = b using the
5347
preconditioned bi-conjugate gradient stabilized method
5349
The return value indicates convergence within max_iter (input)
5350
iterations (0), or no convergence within max_iter iterations (1).
5351
Upon successful return, output arguments have the following values:
5353
approximate solution to Ax = b
5356
the number of iterations performed before the tolerance was reached
5359
the residual after the final iteration
5364
`bicgstab' is an iterative template routine.
5366
`bicgstab' follows the algorithm described on p. 24 in
5368
Templates for the Solution of Linear Systems: Building
5369
Blocks for Iterative Methods, 2nd Edition, R.
5370
Barrett, M. Berry, T. F. Chan, J. Demmel, J. Donato, J. Dongarra,
5371
V. Eijkhout, R. Pozo, C. Romine, H. Van der Vorst,
5372
SIAM, 1994, `ftp.netlib.org/templates/templates.ps'.
5374
The present implementation is inspired from `IML++ 1.2' iterative
5375
method library, `http://math.nist.gov/iml++'.
5378
File: rheolef.info, Node: puzawa algorithm, Up: Algorithms
5380
6.9 `puzawa' - Uzawa algorithm.
5381
===============================
5383
(Source file: `skit/lib/puzawa.h')
5388
template <class Matrix, class Vector, class Preconditioner, class Real>
5389
int puzawa (const Matrix &A, Vector &x, const Vector &b, const Preconditioner &M,
5390
int &max_iter, Real &tol, const Real& rho, std::ostream *p_cerr=0);
5395
The simplest call to 'puzawa' has the folling form:
5396
size_t max_iter = 100;
5398
int status = puzawa(A, x, b, EYE, max_iter, tol, 1.0, &cerr);
5403
`puzawa' solves the linear system A*x=b using the Uzawa method. The
5404
Uzawa method is a descent method in the direction opposite to the
5405
gradient, with a constant step length 'rho'. The convergence is
5406
assured when the step length 'rho' is small enough. If matrix A is
5407
symmetric positive definite, please uses 'pcg' that computes
5408
automatically the optimal descdnt step length at each iteration.
5410
The return value indicates convergence within max_iter (input)
5411
iterations (0), or no convergence within max_iter iterations (1).
5412
Upon successful return, output arguments have the following values:
5414
approximate solution to Ax = b
5417
the number of iterations performed before the tolerance was reached
5420
the residual after the final iteration
5425
template < class Matrix, class Vector, class Preconditioner, class Real, class Size>
5426
int puzawa(const Matrix &A, Vector &x, const Vector &Mb, const Preconditioner &M,
5427
Size &max_iter, Real &tol, const Real& rho,
5428
std::ostream *p_cerr, std::string label)
5430
Vector b = M.solve(Mb);
5431
Real norm2_b = dot(Mb,b);
5432
Real norm2_r = norm2_b;
5433
if (norm2_b == Real(0)) norm2_b = 1;
5434
if (p_cerr) (*p_cerr) << "[" << label << "] #iteration residue" << std::endl;
5435
for (Size n = 0; n <= max_iter; n++) {
5436
Vector Mr = A*x - Mb;
5437
Vector r = M.solve(Mr);
5438
norm2_r = dot(Mr, r);
5439
if (p_cerr) (*p_cerr) << "[" << label << "] " << n << " " << sqrt(norm2_r/norm2_b) << std::endl;
5440
if (norm2_r <= sqr(tol)*norm2_b) {
5441
tol = sqrt(norm2_r/norm2_b);
5447
tol = sqrt(norm2_r/norm2_b);
5452
File: rheolef.info, Node: gmres algorithm, Up: Algorithms
5454
6.10 `gmres' - generalized minimum residual method
5455
==================================================
5457
(Source file: `skit/lib/gmres.h')
5462
template <class Operator, class Vector, class Preconditioner,
5463
class Matrix, class Real, class Int>
5464
int gmres (const Operator &A, Vector &x, const Vector &b,
5465
const Preconditioner &M, Matrix &H, Int m, Int &max_iter, Real &tol);
5470
The simplest call to `gmres' has the folling form:
5473
int status = gmres(a, x, b, EYE, H, m, 100, 1e-7);
5478
`gmres' solves the unsymmetric linear system Ax = b using the
5479
generalized minimum residual method.
5481
The return value indicates convergence within max_iter (input)
5482
iterations (0), or no convergence within max_iter iterations (1).
5483
Upon successful return, output arguments have the following values:
5485
approximate solution to Ax = b
5488
the number of iterations performed before the tolerance was reached
5491
the residual after the final iteration
5492
In addition, M specifies a preconditioner, H specifies a matrix to
5493
hold the coefficients of the upper Hessenberg matrix constructed by
5494
the `gmres' iterations, `m' specifies the number of iterations for
5497
`gmres' requires two matrices as input, A and H. The matrix A, which
5498
will typically be a sparse matrix) corresponds to the matrix in the
5499
linear system Ax=b. The matrix H, which will be typically a dense
5500
matrix, corresponds to the upper Hessenberg matrix H that is
5501
constructed during the `gmres' iterations. Within `gmres', H is used
5502
in a different way than A, so its class must supply different
5503
functionality. That is, A is only accessed though its matrix-vector
5504
and transpose-matrix-vector multiplication functions. On the other
5505
hand, `gmres' solves a dense upper triangular linear system of
5506
equations on H. Therefore, the class to which H belongs must provide
5507
H(i,j) operator for element acess.
5512
It is important to remember that we use the convention that indices
5513
are 0-based. That is H(0,0) is the first component of the matrix H.
5514
Also, the type of the matrix must be compatible with the type of
5515
single vector entry. That is, operations such as H(i,j)*x(j) must be
5516
able to be carried out.
5518
`gmres' is an iterative template routine.
5520
`gmres' follows the algorithm described on p. 20 in
5522
Templates for the Solution of Linear Systems: Building
5523
Blocks for Iterative Methods, 2nd Edition, R.
5524
Barrett, M. Berry, T. F. Chan, J. Demmel, J. Donato, J. Dongarra,
5525
V. Eijkhout, R. Pozo, C. Romine, H. Van der Vorst,
5526
SIAM, 1994, `ftp.netlib.org/templates/templates.ps'.
5528
The present implementation is inspired from `IML++ 1.2' iterative
5529
method library, `http://math.nist.gov/iml++'.
5534
template < class Matrix, class Vector, class Int >
5536
Update(Vector &x, Int k, Matrix &h, Vector &s, Vector v[])
5541
for (Int i = k; i >= 0; i--) {
5543
for (Int j = i - 1; j >= 0; j--)
5544
y(j) -= h(j,i) * y(i);
5547
for (Int j = 0; j <= k; j++)
5552
template < class Real >
5556
return (x > Real(0) ? x : -x);
5561
template < class Operator, class Vector, class Preconditioner,
5562
class Matrix, class Real, class Int >
5564
gmres(const Operator &A, Vector &x, const Vector &b,
5565
const Preconditioner &M, Matrix &H, const Int &m, Int &max_iter,
5570
Vector s(m+1), cs(m+1), sn(m+1), w;
5572
Real normb = norm(M.solve(b));
5573
Vector r = M.solve(b - A * x);
5574
Real beta = norm(r);
5576
if (normb == Real(0))
5579
if ((resid = norm(r) / normb) <= tol) {
5585
Vector *v = new Vector[m+1];
5587
while (j <= max_iter) {
5588
v[0] = r * (1.0 / beta); // ??? r / beta
5592
for (i = 0; i < m && j <= max_iter; i++, j++) {
5593
w = M.solve(A * v[i]);
5594
for (k = 0; k <= i; k++) {
5595
H(k, i) = dot(w, v[k]);
5596
w -= H(k, i) * v[k];
5598
H(i+1, i) = norm(w);
5599
v[i+1] = w * (1.0 / H(i+1, i)); // ??? w / H(i+1, i)
5601
for (k = 0; k < i; k++)
5602
ApplyPlaneRotation(H(k,i), H(k+1,i), cs(k), sn(k));
5604
GeneratePlaneRotation(H(i,i), H(i+1,i), cs(i), sn(i));
5605
ApplyPlaneRotation(H(i,i), H(i+1,i), cs(i), sn(i));
5606
ApplyPlaneRotation(s(i), s(i+1), cs(i), sn(i));
5608
if ((resid = abs(s(i+1)) / normb) < tol) {
5609
Update(x, i, H, s, v);
5616
Update(x, m - 1, H, s, v);
5617
r = M.solve(b - A * x);
5619
if ((resid = beta / normb) < tol) {
5631
template<class Real>
5632
void GeneratePlaneRotation(Real &dx, Real &dy, Real &cs, Real &sn)
5634
if (dy == Real(0)) {
5637
} else if (abs(dy) > abs(dx)) {
5638
Real temp = dx / dy;
5639
sn = 1.0 / ::sqrt( 1.0 + temp*temp );
5642
Real temp = dy / dx;
5643
cs = 1.0 / ::sqrt( 1.0 + temp*temp );
5647
template<class Real>
5648
void ApplyPlaneRotation(Real &dx, Real &dy, Real &cs, Real &sn)
5650
Real temp = cs * dx + sn * dy;
5651
dy = -sn * dx + cs * dy;
5656
File: rheolef.info, Node: Forms, Up: Top
5679
File: rheolef.info, Node: convect form, Up: Forms
5681
7.1 `convect' - discontinuous Galerkin
5682
======================================
5684
(Source file: `nfem/form_element/convect.h')
5689
form(const space& V, const space& V, "convect", uh);
5694
Assembly the matrix associated to the discontinuous Galerkin method
5695
on the finite element space V.
5698
c(u,v) = | a.grad(u) v dx + skips...
5701
for all u,v in V, where the vector field a is given. The V
5702
space may `P1d' finite element spaces for building the form.
5705
File: rheolef.info, Node: d_dx form, Up: Forms
5707
7.2 `d_dx'I - derivatives
5708
=========================
5710
(Source file: `nfem/form_element/d_dx0.h')
5715
form (const space V, const space& M, "d_dx0");
5716
form (const space V, const space& M, "d_dx1");
5717
form (const space V, const space& M, "d_dx2");
5722
Assembly the form associated to a derivative operator from the `V'
5723
finite element space to the `M' one:
5726
b_i(u,q) = | ---- q dx, i = 0,1,2
5729
In the axisymetric `rz' case, the form is defined by
5732
b_0(u,q) = | --- q r dr dz
5735
If the V space is a `P1' (resp. `P2') finite element space, the M
5736
space may be a `P0' (resp. `P1d') one.
5741
The following piece of code build the Laplacian form associated to
5742
the P1 approximation:
5743
geo omega ("square");
5744
space V (omega, "P1");
5745
space M (omega, "P0");
5746
form b (V, M, "d_dx0");
5751
Only edge, triangular and tetrahedal finite element meshes are yet
5755
File: rheolef.info, Node: 2W form, Up: Forms
5757
7.3 `2W' - vorticity tensor
5758
===========================
5760
(Source file: `nfem/form_element/2W.h')
5765
form (const space V, const space& T, "2W");
5770
Assembly the form associated to the vorticity tensor, i.e. the
5771
unsymmetric part of the gradient of a vector field. These
5772
derivative are usefull in fluid mechanic.
5775
b(u,tau) = | 2 W(u) : tau dx,
5779
2 D(u) = grad u - (grad u)^T
5780
If the V space is a vector-valued `P1' (resp. `P2') finite element
5781
space, the T space may be a tensor-valued `P0' (resp. `P1d') one.
5786
The following piece of code build the Laplacian form associated to
5787
the P1 approximation:
5788
geo omega ("square");
5789
space V (omega, "P1", "vector");
5790
space T (omega, "P0", "tensor");
5791
form b (V, T, "2W");
5794
File: rheolef.info, Node: curl form, Up: Forms
5796
7.4 `curl' - curl operator
5797
==========================
5799
(Source file: `nfem/form_element/curl.h')
5804
form(const space V, const space& M, "curl");
5809
Assembly the form associated to the curl operator on finite element
5810
space. In three dimensions, both V and M are vector-valued:
5813
b(u,q) = | curl(u).q dx
5816
In two dimensions, only V is vector-valued. The V space may be
5817
a either `P2' finite element space, while the M space may be
5818
`P1d'. See also *note form class:: and *note space class::.
5823
The following piece of code build the divergence form associated to
5824
the `P2' approximation for a three dimensional geometry:
5826
space V(omega, "P2", "vector");
5827
space M(omega, "P1d", "vector");
5828
form b(V, M, "curl");
5829
while this code becomes in two dimension:
5830
geo omega("square");
5831
space V(omega, "P2", "vector");
5832
space M(omega, "P1d");
5833
form b(V, M, "curl");
5836
File: rheolef.info, Node: 2D_D form, Up: Forms
5838
7.5 `2D_D' - -div 2D operator
5839
=============================
5841
(Source file: `nfem/form_element/2D_D.h')
5846
form (const space V, const space& V, "2D_D");
5851
Assembly the form associated to the `div(2D(.))' components of the
5852
operator on the finite element space `V':
5855
a(u,v) = | 2 Dij(ui) Dij(vj) dx
5860
Dij(u) = - ( ---- + ---- )
5862
This form is usefull when considering elasticity or Stokes problems.
5867
Here is an example of the vector-valued form:
5868
geo omega ("square");
5869
space V (omega, "P2", "vector");
5870
form a (V, V, "2D_D");
5871
Note that a factor two is here applied to the form. This factor is
5872
commonly used in practice.
5875
File: rheolef.info, Node: d_ds form, Up: Forms
5877
7.6 `d_ds' - Curvilinear derivative
5878
===================================
5880
(Source file: `nfem/form_element/d_ds.h')
5885
form(const space& V, const space& W, "d_ds");
5890
Assembly the matrix associated to the following integral:
5893
m(u,v) = | du/ds v ds
5898
File: rheolef.info, Node: d2_ds2 form, Up: Forms
5900
7.7 `d2_ds2' - Curvilinear second derivative
5901
============================================
5903
(Source file: `nfem/form_element/d2_ds2.h')
5908
form(const space& V, const space& W, "d2_ds2");
5913
Assembly the matrix associated to the following integral:
5916
m(u,v) = | d^2u/ds^2 v ds
5921
File: rheolef.info, Node: div form, Up: Forms
5923
7.8 `div' - divergence operator
5924
===============================
5926
(Source file: `nfem/form_element/div.h')
5931
form(const space V, const space& M, "div");
5936
Assembly the form associated to the divergence operator on a finite
5940
b(u,q) = | div(u) q dx
5943
The V space may be a either `P1' or `P2' finite element space,
5944
while the M space may be `P0' or `P1d' respectively. See also
5945
*note form class:: and *note space class::.
5950
The following piece of code build the divergence form associated to
5951
the `P1' approximation:
5952
geo omega("square");
5953
space V(omega, "P1", "vector");
5954
space M(omega, "P0");
5955
form b(V, M, "div");
5958
File: rheolef.info, Node: mass form, Up: Forms
5960
7.9 `mass' - L2 scalar product
5961
==============================
5963
(Source file: `nfem/form_element/mass.h')
5968
form(const space& V, const space& V, "mass");
5969
form(const space& M, const space& V, "mass");
5970
form (const space& V, const space& V, "mass", const domain& gamma);
5971
form_diag(const space& V, "mass");
5976
Assembly the matrix associated to the L2 scalar product of the
5977
finite element space V.
5984
The V space may be either a `P0', `P1', `P2', `bubble', `P1d'
5985
and `P1d' finite element spaces for building a form
5988
The use of quadrature formulae is sometime usefull for building
5989
diagonal matrix. These approximate matrix are eay to invert.
5990
This procedure is available for `P0' and `P1' approximations.
5992
Notes that when dealing with discontinuous finite element space,
5993
i.e. `P0' and `P1d', the corresponding mass matrix is block
5994
diagonal, and the `inv_mass' form may be usefull.
5996
When two different space M and V are supplied, assembly the matrix
5997
associated to the projection operator from one finite element
6004
for all q in M and v in V.
6006
This form is usefull for instance to convert discontinuous gradient
6007
components to a continuous approximation. The transpose operator may
6008
also be usefull to performs the opposite operation.
6010
The following $V$ and $M$ space approximation combinations are
6011
supported for the `mass' form: P0-P1, P0-P1d, P1d-P2,
6017
The following piece of code build the mass matrix associated to the
6021
form m(V, V, "mass");
6022
The use of lumped mass form write also:
6023
form_diag md(V, "mass");
6024
The following piece of code build the projection form:
6028
form m(M, V, "mass");
6030
Scalar product on the boundary
6031
------------------------------
6033
Assembly the matrix associated to the L2 scalar product related to
6034
a boundary domain of a mesh and a specified polynomial
6035
approximation. These forms are usefull when defining
6036
non-homogeneous Neumann or Robin boundary conditions.
6038
Let W be a space of functions defined on Gamma, a subset of the
6039
boundary of the whole domain Omega.
6045
for all u, v in W. Let V a space of functions defined on Omega
6046
and gamma the trace operator from V into W. For all u in W and
6050
mb(u,v) = | u gamma(v) dx
6053
For all u and v in V:
6056
ab(u,v) = | gamma(u) gamma(v) dx
6063
The following piece of code build forms for the P1 approximation,
6064
assuming that the mesh contains a domain named `boundary':
6065
geo omega ("square");
6066
domain gamma = omega.boundary();
6067
space V (omega, "P1");
6068
space W (omega, gamma, "P1");
6069
form m (W, W, "mass");
6070
form mb (W, V, "mass");
6071
form ab (V, V, "mass", gamma);
6074
File: rheolef.info, Node: grad form, Up: Forms
6076
7.10 `grad' - gradient operator
6077
===============================
6079
(Source file: `nfem/form_element/grad.h')
6084
form(const space V, const space& M, "grad");
6089
Assembly the form associated to the gradient operator on a finite
6093
b(u, q) = | grad(u).q dx
6096
The V space may be a either `P1' or `P2' finite element space,
6097
while the M space may be `P0' or `P1d' respectively. See also
6098
*note form class:: and *note space class::.
6103
The following piece of code build the divergence form associated to
6104
the `P1' approximation:
6105
geo omega("square");
6106
space V(omega, "P1");
6107
space M(omega, "P0", "vector");
6108
form b(V, M, "grad");
6111
File: rheolef.info, Node: inv_mass form, Up: Forms
6113
7.11 `inv_mass' - invert of L2 scalar product
6114
=============================================
6116
(Source file: `nfem/form_element/inv_mass.h')
6121
form(const space& V, const space& V, "inv_mass");
6126
Assembly the invert of the matrix associated to the L2 scalar product
6127
of the finite element space V:
6133
The V space may be either a `P0' or `P1d' discontinuous finite
6134
element spaces *note form class::.
6139
The following piece of code build the invert of the mass matrix
6140
associated to the `P1d' approximation:
6141
geo omega_h ("square");
6142
space Vh (omega_h , "P1d");
6143
form im (Vh, Vh, "inv_mass");
6146
File: rheolef.info, Node: grad_grad form, Up: Forms
6148
7.12 `grad_grad' - -Laplacian operator
6149
======================================
6151
(Source file: `nfem/form_element/grad_grad.h')
6156
form(const space V, const space& V, "grad_grad");
6161
Assembly the form associated to the Laplacian operator on a finite
6165
a(u,v) = | grad(u).grad(v) dx
6168
The V space may be a either `P1' , `P2' or `P1d' finite
6169
element space. See also *note form class:: and *note space class::.
6174
The following piece of code build the Laplacian form associated to
6175
the `P1' approximation:
6178
form a(V, V, "grad_grad");
6181
File: rheolef.info, Node: div_div form, Up: Forms
6183
7.13 `div_div' - -grad div operator
6184
===================================
6186
(Source file: `nfem/form_element/div_div.h')
6191
form (const space V, const space& V, "div_div");
6196
Assembly the form associated to the `grad(div(.))' operator on the
6197
finite element space `V':
6200
a(u,v) = | div(u) div(v) dx
6203
This form is usefull when considering elasticity problem.
6208
Here is an example of the vector-valued form:
6209
geo omega ("square");
6210
space V (omega, "P2", "vector");
6211
form a (V, V, "div_div");
6214
File: rheolef.info, Node: 2D form, Up: Forms
6216
7.14 `2D' - rate of deformation tensor
6217
======================================
6219
(Source file: `nfem/form_element/2D.h')
6224
form (const space V, const space& T, "2D");
6229
Assembly the form associated to the rate of deformation tensor, i.e.
6230
the symmetric part of the gradient of a vector field. These
6231
derivative are usefull in fluid mechanic and elasticity.
6234
b(u,tau) = | 2 D(u) : tau dx,
6238
2 D(u) = grad u + (grad u)^T
6240
If the V space is a vector-valued `P1' (resp. `P2') finite element
6241
space, the T space may be a tensor-valued `P0' (resp. `P1d') one.
6246
The following piece of code build the Laplacian form associated to
6247
the P1 approximation:
6248
geo omega ("square");
6249
space V (omega, "P1", "vector");
6250
space T (omega, "P0", "tensor");
6251
form b (V, T, "2D");
6254
File: rheolef.info, Node: Internals, Up: Top
6262
* form_element internal::
6264
* characteristic internal::
6267
* reference_element internal::
6268
* quadrature internal::
6269
* quadrangle internal::
6272
* triangle internal::
6273
* geo_element internal::
6276
* numbering internal::
6277
* smart_pointer internal::
6278
* heap_allocator internal::
6279
* stack_allocator internal::
6280
* pretty_name internal::
6281
* acinclude internal::
6284
File: rheolef.info, Node: geomap internal, Up: Internals
6286
8.1 `geomap' - discrete mesh advection by a field: gh(x)=fh(x-dt*uh
6287
===================================================================
6289
(Source file: `nfem/lib/geomap.h')
6294
The class geomap is a fundamental class used for the correspondance
6295
between fields defined on different meshes or for advection problems.
6296
This class is used for the method of characteristic.
6301
The following code compute gh(x)=fh(x-dt*uh(x)):
6302
field uh = interpolate (Vh, u);
6303
field fh = interpolate (Fh, f);
6304
geomap X (Fh, -dt*uh);
6305
field gh = compose (fh, X);
6306
For a complete example, see `convect.cc' in the example directory.
6311
struct geomap_option_type {
6312
size_t n_track_step; // loop while tracking: y = X_u(x)
6313
geomap_option_type() : n_track_step(1) {}
6315
class geomap : public Vector<meshpoint> {
6317
geomap () : advected(false), use_space(true) {}
6318
// Maps a quadrature-dependent lattice over Th_1 onto Th_2 triangulation.
6319
geomap (const geo& Th_2, const geo& Th_1,
6320
std::string quadrature="default", size_t order=0, bool allow_approximate_edges=true);
6322
// Maps a quadrature-dependent lattice over Th_1 onto Th_2 triangulation
6323
// thro' an offset vector field u to be defined on Th_1.
6324
geomap (const geo& Th_2, const geo& Th_1, const field& advection_h_1,
6325
std::string quadrature="default", size_t order=0) ;
6327
// TO BE REMOVED: backward compatibility
6328
// Maps space Vh_1 dof's onto the meshpoints of Th_2 triangulation
6329
/* geomap : ( Vh_1.dof's ) |--> ( Th_2 )
6331
geomap (const geo& Th_2, const space& Vh_1,
6332
bool allow_approximate_edges=true );
6333
// Same but for P1d, using points inside the triangle to preserve discontinuity
6334
geomap (const geo& Th_2, const space& Vh_1,
6335
Float barycentric_weight );
6337
// TO BE REMOVED: backward compatibility
6338
// Maps space Vh_1 dof's onto the meshpoints of Th_2 triangulation
6339
// thro' an offset u to be defined on Vh_1 dof's.
6340
/* geomap : ( Vh_1.dof's ) |--> ( Th_2 )
6342
geomap (const geo& Th_2, const space& Vh_1, const field& advection_h_1);
6344
// Does the same, but assumes that Th_2 is the triangulation for Vh_1
6345
geomap (const space& Vh_2, const field& advection_h_1, const geomap_option_type& opts = geomap_option_type());
6350
// TO BE REMOVED: backward compatibility
6353
{ if (!use_space) error_macro("Lattice-based geomaps use no space");
6358
get_origin_triangulation () const
6362
get_target_triangulation () const
6368
// TO BE REMOVED: backward compatibility
6376
// TO BE REMOVED: backward compatibility
6378
return _quadrature; }
6380
bool is_inside (size_t dof) const { return _is_inside [dof]; }
6381
bool no_barycentric_weight() const { return _barycentric_weight == Float(0); }
6382
Float barycentric_weight() const { return _barycentric_weight; }
6385
friend class fieldog;
6387
void init (const field& advection);
6388
// non use_space mode
6390
meshpoint advect (const point& q, size_t iK_Th_1);
6391
meshpoint robust_advect_1 (const point& x0, const point& va, bool& is_inside) const;
6392
meshpoint robust_advect_N (const point& x0, const point& v0, const field& vh, size_t n, bool& is_inside) const;
6400
std::string _quadrature;
6402
bool _allow_approximate_edges;
6408
std::vector<point> quad_point;
6409
std::vector<Float> quad_weight;
6411
// flag when a dof go outside of the domain
6412
std::vector<bool> _is_inside;
6414
Float _barycentric_weight;
6415
geomap_option_type _option;
6419
geomap::geomap (const geo& Th_2, const space& Vh_1, const field& advection) :
6420
_Th_1 (Vh_1.get_geo()), _Th_2 (Th_2), _advection(advection), advected(true),
6421
_quadrature("default"), _order(0),
6422
_allow_approximate_edges(true), use_space(true), _Vh_1(Vh_1),
6423
_is_inside(), _barycentric_weight(0)
6424
{ init (advection); };
6427
geomap::geomap (const space& Vh_1, const field& advection, const geomap_option_type& opt) :
6428
_Th_1 (Vh_1.get_geo()), _Th_2 (Vh_1.get_geo()), _advection(advection), advected(true),
6429
_quadrature("default"), _order(0),
6430
_allow_approximate_edges(true), use_space(true), _Vh_1(Vh_1),
6431
_is_inside(), _barycentric_weight(0), _option(opt)
6432
{ init (advection); };
6438
const field& advection,
6439
std::string quadrature,
6444
_advection(advection),
6446
_quadrature(quadrature),
6448
_allow_approximate_edges(true),
6452
_barycentric_weight(0),
6455
if (_advection.get_geo() !=_Th_1)
6456
error_macro("The advection field should be defined on the original mesh Th_1="
6457
<< Th_1.name() << " and not " << _advection.get_geo().name());
6465
std::string quadrature,
6467
bool allow_approximate_edges)
6472
_quadrature(quadrature),
6474
_allow_approximate_edges(allow_approximate_edges),
6478
_barycentric_weight(0),
6484
// Composition with a field
6485
/* f being a field of space V, X a geomap on this space, foX=compose(f,X) is their
6488
field compose (const field& f, const geomap& g);
6489
/* send also a callback function when advection goes outside a domain
6490
* this is usefull when using upstream boundary conditions
6492
field compose (const field& f, const geomap& g, Float (*f_outside)(const point&));
6493
template <class Func>
6494
field compose (const field& f, const geomap& g, Func f_outside);
6499
class fieldog: public field // and friend of geomap.
6503
fieldog (const field& _f, const geomap& _g) : f(_f), g(_g) {};
6507
operator() (const point& x) const;
6510
operator() (const meshpoint& x) const;
6513
operator() (const point& x_hat, size_t e) const
6514
{ return operator() (meshpoint(x_hat,e)); };
6517
quadrature_values (size_t iK) const;
6520
quadrature_d_dxi_values (size_t i, size_t iK) const;
6524
{ return g._quadrature; };
6528
{ return g._order; };
6532
{ return f.get_space(); };
6540
File: rheolef.info, Node: form_element internal, Up: Internals
6542
8.2 `form_element' - bilinear form on a single element
6543
======================================================
6545
(Source file: `nfem/lib/form_element.h')
6550
The `form_element' class defines functions that compute a bilinear
6551
form defined between two polynomial basis on a single geometrical
6552
element. This bilinear form is represented by a matrix.
6554
The bilinear form is designated by a string, e.g. "mass", "grad_grad",
6555
... indicating the form. The form depends also of the geometrical
6556
element: triangle, square, tetrahedron (*note geo_element internal::).
6561
The `form_element' class is managed by (*note smart_pointer
6562
internal::). This class uses a pointer on a pure virtual class
6563
`form_element_rep' while the effective code refers to the specific
6564
concrete derived classes: mass, grad_grad, etc.
6569
class form_element : public smart_pointer<form_element_rep> {
6574
typedef size_t size_type;
6578
form_element (std::string name = "");
6582
void initialize (const space& X, const space& Y) const;
6584
// optional, for scalar-weighted forms:
6585
void set_weight (const field& wh) const;
6587
// for special forms for which coordinate-system specific weights (ie, cylindrical
6588
// coordinates weights) should not be used:
6589
void set_use_coordinate_system_weight(bool use) const;
6593
void operator() (const geo_element& K, ublas::matrix<Float>& m) const;
6594
const field& get_weight() const;
6598
File: rheolef.info, Node: iofem internal, Up: Internals
6600
8.3 `iofem' - input and output finite element manipulators
6601
==========================================================
6603
(Source file: `nfem/lib/iofem.h')
6608
This class implements some specific finite element manipulators.
6609
For a general presentation of stream manipulators, reports to class
6612
Valuated manipulators
6613
---------------------
6617
set a cutting plane for visualizations and post-processing.
6618
cout << cut << origin(point(0.5, 0.5, 0.5)) << normal(point(1,1,0) << uh;
6622
specifies the topography field when representing a bidimensionnal
6623
field in tridimensionnal elevation.
6624
cout << topography(zh) << uh;
6625
This manipulator takes an agument that specifies a scalar
6626
field value `zh' for the elevation, while `uh' is the
6627
scalar field to represent. Then, the z-elevation takes
6631
File: rheolef.info, Node: characteristic internal, Up: Internals
6633
8.4 `characteristic' - discrete mesh advection by a field: gh(x)=fh(x-dt*uh
6634
===========================================================================
6636
(Source file: `nfem/lib/characteristic.h')
6641
The class characteristic is a class implementing the Lagrange-Galerkin
6642
method. It is the extension of the method of characteristic in the
6643
finite element context. This is an expreimental implementation:
6644
please use the "geomap" one for practical usage.
6649
The following code compute the Riesz representant, denoted by "mgh" of
6650
gh(x)=fh(x-dt*uh(x)).
6654
characteristic X (omega_h, -dt*uh);
6655
field mgh = riesz_representer(Vh, compose(fh, X));
6656
The Riesz representer is the "mgh" vector of values:
6657
mgh(i) = integrate fh(x-dt*uh(x)) phi_i(x) dx
6658
where phi_i is the i-th basis function in Vh and the integral is
6659
evaluated by using a quadrature formulae. By default the quadrature
6660
formule is Gauss-Lobatto with the order equal to the polynomial order
6661
of Vh. This quadrature formulae guaranties inconditional stability
6662
at any polynomial order (order 1: trapeze, order 2: simpson).
6663
Extension will accept in the future alternative quadrature formulae.
6668
class characteristic {
6670
characteristic(const geo& bg_omega, const field& ah);
6671
const geo& get_background_geo() const;
6672
const field& get_advection() const;
6677
class field_o_characteristic {
6679
field_o_characteristic(const field& fh, const characteristic& X);
6680
friend field_o_characteristic compose(
6681
const field& fh, const characteristic& X);
6682
Float operator() (const point& x) const;
6683
point vector_evaluate (const point& x) const;
6684
tensor tensor_evaluate (const point& x) const;
6685
const field& get_field() const;
6686
const geo& get_background_geo() const;
6687
const field& get_advection() const;
6691
point advect (const point& x0) const;
6695
File: rheolef.info, Node: point internal, Up: Internals
6697
8.5 `point' - Point reference element
6698
=====================================
6700
(Source file: `nfem/basis/point.icc')
6705
The point reference element is defined for convenience. It is a
6706
0-dimensional element with measure equal to 1.
6711
const size_t dimension = 0;
6712
const size_t measure = 1;
6715
File: rheolef.info, Node: basis internal, Up: Internals
6717
8.6 `basis' - polynomial basis
6718
==============================
6720
(Source file: `nfem/basis/basis.h')
6725
The `basis' class defines functions that evaluates a polynomial basis
6726
and its derivatives on a point. The polynomial basis is designated by
6727
a string, e.g. "P0", "P1", "P2", "bubble",... indicating the basis.
6728
The basis depends also of the reference element: triangle, square,
6729
tetrahedron (*note reference_element internal::). For instance, on a
6730
square, the "P1" string designates the common Q1 four-nodes basis on
6731
the reference square.
6733
The nodes associated to the Lagrange polynomial basis are also
6734
available by its associated accessor.
6739
The `basis' class is a *note smart_pointer internal::) class on a
6740
`basis_rep' class that is a pure virtual base class for effective
6741
bases, e.g. basis_P1, basis_P1, etc.
6746
class basis : public smart_pointer<basis_rep> {
6751
typedef size_t size_type;
6752
typedef basis_rep::dof_family_type dof_family_type;
6756
basis (std::string name = "");
6760
std::string name() const;
6761
size_type degree() const;
6762
size_type size (reference_element hat_K, dof_family_type family=reference_element::dof_family_max) const;
6763
dof_family_type family() const;
6765
dof_family_type dof_family(
6766
reference_element hat_K,
6767
size_type i_dof_local) const;
6770
reference_element hat_K,
6771
size_type i_dof_local,
6772
const point& hat_x) const;
6775
reference_element hat_K,
6776
size_type i_dof_local,
6777
const point& hat_x) const;
6779
basic_point<point> hessian_eval(
6780
reference_element hat_K,
6781
size_type i_dof_local,
6782
const point& hat_x) const;
6785
reference_element hat_K,
6787
std::vector<Float>& values) const;
6790
reference_element hat_K,
6792
std::vector<point>& values) const;
6795
reference_element hat_K,
6797
std::vector<basic_point<point> >& values) const;
6800
reference_element hat_K,
6801
std::vector<point>& hat_node) const;
6803
void dump(std::ostream& out = std::cerr) const;
6807
File: rheolef.info, Node: reference_element internal, Up: Internals
6809
8.7 `reference_element' - reference element
6810
===========================================
6812
(Source file: `nfem/basis/reference_element.h')
6817
The `reference_element' class defines all supported types of
6818
geometrical elements in one, two and three dimensions. The set of
6819
supported elements are designate by a letter
6827
triangle(dimension 2)
6830
quadrangle(dimension 2)
6833
tetrahedron(dimension 3)
6839
hexaedron(dimension 3)
6844
class reference_element {
6850
typedef std::vector<int>::size_type size_type;
6852
// defines enum_type { p, t, q ..., H, ...};
6853
// in an automatically generated file :
6874
static const size_type not_set;
6875
static const size_type max_n_subgeo = 12;
6876
static const size_type max_subgeo_vertex = 8;
6878
// allocators/deallocators:
6880
reference_element (enum_type x = max_size);
6884
enum_type type() const;
6886
size_type dimension() const;
6887
friend Float measure (reference_element hat_K);
6888
size_type size() const;
6889
size_type n_edge() const;
6890
size_type n_face() const;
6891
size_type n_subgeo(size_type subgeo_dim) const;
6892
size_type subgeo_size(size_type subgeo_dim, size_type i_subgeo) const;
6893
size_type subgeo_local_vertex(size_type subgeo_dim, size_type i_subgeo, size_type i_subgeo_vertex) const;
6894
size_type heap_size() const;
6895
size_type heap_offset(size_type subgeo_dim) const;
6897
void set_type (enum_type x);
6898
void set_type (size_type n_vertex, size_type dim);
6899
void set_name (char name);
6908
// these constants are initialized in
6909
// "reference_element_declare.c"
6910
// automatically generated.
6912
static const char _name [max_size];
6913
static const size_type _dimension [max_size];
6914
static const size_type _n_subgeo [max_size][4];
6915
static const size_type _subgeo_size [max_size*4*max_n_subgeo];
6916
static const size_type _subgeo_local_vertex [max_size*4*max_n_subgeo*max_subgeo_vertex];
6917
static const size_type _heap_size [max_size];
6918
static const size_type _heap_offset [max_size][4];
6919
friend std::istream& operator>>(std::istream&, class geo_element&);
6923
File: rheolef.info, Node: quadrature internal, Up: Internals
6925
8.8 `quadrature' - quadrature formulae on the reference lement
6926
==============================================================
6928
(Source file: `nfem/basis/quadrature.h')
6933
The `quadrature' class defines a container for a quadrature formulae on
6934
the reference element (*note reference_element internal::). This
6935
container stores the nodes coordinates and the weights.
6937
The constructor takes two arguments
6938
-----------------------------------
6940
the reference element K and the order r of the quadrature formulae.
6941
The formulae is exact when computing the integral of a polynom p that
6942
degree is less or equal to order r.
6945
| p(x) dx = \ p(x_q) w_q
6952
The formulae is optimal when it uses a minimal number of nodes n.
6953
Optimal quadrature formula are hard-coded in this class. Not all
6954
reference elements and orders are yet implemented. This class will be
6955
completed in the future.
6965
typedef quadrature_on_geo::size_type size_type;
6966
typedef quadrature_option_type::family_type family_type;
6967
typedef std::vector<weighted_point>::const_iterator const_iterator;
6971
quadrature (quadrature_option_type opt = quadrature_option_type());
6975
void set_order (size_type order);
6976
void set_family (family_type ft);
6980
size_type get_order() const;
6981
family_type get_family() const;
6982
std::string get_family_name() const;
6983
size_type size (reference_element hat_K) const;
6984
const_iterator begin (reference_element hat_K) const;
6985
const_iterator end (reference_element hat_K) const;
6986
friend std::ostream& operator<< (std::ostream&, const quadrature&);
6988
quadrature_option_type _options;
6989
mutable quadrature_on_geo _quad [reference_element::max_size];
6990
mutable std::vector<bool> _initialized;
6991
void _initialize (reference_element hat_K) const;
6993
quadrature (const quadrature&);
6994
quadrature operator= (const quadrature&);
6998
File: rheolef.info, Node: quadrangle internal, Up: Internals
7000
8.9 `quadrangle' - Quadrangular reference element
7001
=================================================
7003
(Source file: `nfem/basis/quadrangle.icc')
7008
The quadrangular reference element is [0,1]^2.
7026
const size_t dimension = 2;
7027
const Float measure = 4;
7028
const size_t n_vertex = 4;
7029
const Float vertex [n_vertex][dimension] = {
7034
const size_t n_edge = 4;
7035
const size_t edge [n_edge][2] = {
7042
File: rheolef.info, Node: prism internal, Up: Internals
7044
8.10 `prism' - Prism reference element
7045
======================================
7047
(Source file: `nfem/basis/prism.icc')
7052
The prism reference element is
7053
K = { 0 < x < 1 and 0 < y < 1-x and -1 < z < 1 }
7058
The orientation is such that triedra (01, 02, 03) is direct and all
7059
faces, see from exterior, are in the direct sens. References: P. L.
7060
Georges, "Generation automatique de maillages", page 24-, coll RMA, 16,
7078
const size_t dimension = 3;
7079
const Float measure = 1;
7080
const size_t n_vertex = 6;
7081
const Float vertex [n_vertex][dimension] = {
7088
const size_t n_face = 5;
7089
const size_t face [n_face][4] = {
7090
{ 0, 2, 1, size_t(-1) },
7093
{ 3, 4, 5, size_t(-1) },
7095
const size_t n_edge = 9;
7096
const size_t edge [n_edge][2] = {
7108
File: rheolef.info, Node: hexa internal, Up: Internals
7110
8.11 `hexa' - Hexaedra reference element
7111
========================================
7113
(Source file: `nfem/basis/hexa.icc')
7118
The hexa reference element is [-1,1]^3.
7123
The orientation is such that triedra (01, 03, 04) is direct and all
7124
faces, see from exterior, are in the direct sens. References: P. L.
7125
Georges, "Generation automatique de maillages", page 24-, coll RMA, 16,
7144
const size_t dimension = 3;
7145
const Float measure = 8;
7146
const size_t n_vertex = 8;
7147
const Float vertex [n_vertex][dimension] = {
7156
const size_t n_face = 6;
7157
const size_t face [n_face][4] = {
7164
const size_t n_edge = 12;
7165
const size_t edge [n_edge][2] = {
7180
File: rheolef.info, Node: triangle internal, Up: Internals
7182
8.12 `triangle' - Triangle reference element
7183
============================================
7185
(Source file: `nfem/basis/triangle.icc')
7190
The triangle reference element is
7191
K = { 0 < x < 1 and 0 < y < 1-x }
7209
const size_t dimension = 2;
7210
const Float measure = 0.5;
7211
const size_t n_vertex = 3;
7212
const Float vertex [n_vertex][dimension] = {
7216
const size_t n_edge = 3;
7217
const size_t edge [n_edge][2] = {
7223
File: rheolef.info, Node: geo_element internal, Up: Internals
7225
8.13 `geo_element' - element of a mesh
7226
======================================
7228
(Source file: `nfem/basis/geo_element_v1.h')
7233
Defines geometrical elements and sides as a set of vertice and edge
7234
indexes. This element is obtained after a Piola transformation
7235
from a reference element (*note reference_element internal::).
7236
Indexes are related to arrays of edges and vertices. These arrays
7237
are included in the description of the mesh. Thus, this class is
7238
related of a given mesh instance (*note geo class::).
7243
This is the test of geo_element:
7246
cout << "n_vertices: " << K.size() << endl
7247
<< "n_edges : " << K.n_edges() << endl
7248
<< "dimension : " << K.dimension() << endl << endl;
7249
for(geo_element::size_type i = 0; i < K.size(); i++)
7251
for(geo_element::size_type i = 0; i < K.n_edges(); i++)
7252
K.set_edge(i, i*10+5) ;
7253
cout << "vertices: local -> global" << endl;
7254
for (geo_element::size_type vloc = 0; vloc < K.size(); vloc++)
7255
cout << vloc << "-> " << K[vloc] << endl;
7257
<< "edges: local -> global" << endl;
7258
for (geo_element::size_type eloc = 0; eloc < K.n_edges(); eloc++) {
7259
geo_element::size_type vloc1 = subgeo_local_vertex(1, eloc, 0);
7260
geo_element::size_type vloc2 = subgeo_local_vertex(1, eloc, 1);
7261
cout << eloc << "-> " << K.edge(eloc) << endl
7262
<< "local_vertex_from_edge(" << eloc
7263
<< ") -> (" << vloc1 << ", " << vloc2 << ")" << endl;
7269
class geo_element : public reference_element {
7272
// allocators/deallocators:
7274
geo_element(enum_type t = max_size);
7275
geo_element(const geo_element&);
7276
explicit geo_element (const class tiny_element&);
7277
void copy (const geo_element&);
7278
void copy (const class tiny_element&);
7279
geo_element& operator = (const geo_element&);
7284
size_type index() const;
7285
size_type operator [] (size_type i) const;
7286
size_type side(size_type i_side) const;
7287
void build_side(size_type i_side, geo_element& S) const;
7290
size_type edge (size_type i_edge) const;
7291
size_type face (size_type i_face) const;
7293
size_type subgeo(const geo_element& S) const;
7294
size_type subgeo (size_type subgeo_dim, size_type i_subgeo) const;
7295
size_type subgeo_vertex (size_type subgeo_dim, size_type i_subgeo,
7296
size_type i_subgeo_vertex) const;
7297
void build_subgeo(size_type subgeo_dim, size_type i_subgeo, geo_element& S) const;
7298
size_type subgeo_local_index(const geo_element& S) const;
7302
void set_type (enum_type t);
7303
void set_type (size_type sz, size_type dim);
7304
void set_name (char name);
7306
void set_index(size_type idx);
7307
size_type& operator [] (size_type i);
7308
void set_side (size_type i_side, size_type idx);
7310
void set_edge (size_type i_edge, size_type idx);
7311
void set_face (size_type i_face, size_type idx);
7312
void set_subgeo (size_type subgeo_dim, size_type i_subgeo,
7314
void set_subgeo(const geo_element& S, size_type idx);
7318
friend std::istream& operator >> (std::istream&, geo_element&);
7319
friend std::ostream& operator << (std::ostream&, const geo_element&);
7320
std::ostream& dump(std::ostream& s = std::cerr) const;
7327
// memory management:
7333
File: rheolef.info, Node: tetra internal, Up: Internals
7335
8.14 `tetra' - Tetraedra reference element
7336
==========================================
7338
(Source file: `nfem/basis/tetra.icc')
7343
The tetraedra reference element is
7344
K = { 0 < x < 1 and 0 < y < 1-x and 0 < z < 1-x-y }
7349
The orientation is such that triedra (01, 02, 03) is direct, and all
7350
faces, see from exterior, are in the direct sens. References: P. L.
7351
Georges, "Generation automatique de maillages", page 24-, coll RMA, 16,
7368
const size_t dimension = 3;
7369
const Float measure = Float(1.)/Float(6.);
7370
const size_t n_vertex = 4;
7371
const Float vertex [n_vertex][dimension] = {
7376
const size_t n_face = 4;
7377
const size_t face [n_face][3] = {
7382
const size_t n_edge = 6;
7383
const size_t edge [n_edge][2] = {
7392
File: rheolef.info, Node: edge internal, Up: Internals
7394
8.15 `edge' - Edge reference element
7395
====================================
7397
(Source file: `nfem/basis/edge.icc')
7402
The edge reference element is K = [0,1].
7408
const size_t dimension = 1;
7409
const size_t measure = 1;
7410
const size_t n_vertex = 2;
7411
const Float vertex [n_vertex][dimension] = {
7416
File: rheolef.info, Node: numbering internal, Up: Internals
7418
8.16 `numbering' - global degree of freedom numbering
7419
=====================================================
7421
(Source file: `nfem/basis/numbering.h')
7426
The `numbering' class defines methods that furnish global numbering
7427
of degrees of freedom. This numbering depends upon the degrees of
7428
polynoms on elements and upon the continuity requirement at
7429
inter-element boundary. For instance the "P1" continuous finite
7430
element approximation has one degree of freedom per vertice of the
7431
mesh, while its discontinuous counterpart has dim(basis) times the
7432
number of elements of the mesh, where dim(basis) is the size of the
7433
local finite element basis.
7438
class numbering : public smart_pointer<numbering_rep> {
7442
typedef numbering_rep rep;
7443
typedef smart_pointer<numbering_rep> base;
7444
typedef size_t size_type;
7448
numbering (std::string name = "");
7449
numbering (numbering_rep* ptr);
7451
virtual ~numbering() {}
7455
std::string name() const;
7458
size_type mesh_map_dimension,
7459
const size_type* mesh_n_geo,
7460
const size_type* mesh_n_element) const;
7464
const size_type* mesh_n_geo,
7465
const size_type* mesh_n_element,
7466
const geo_element& K,
7467
size_type i_dof_local) const;
7471
const size_type* mesh_n_geo,
7472
const size_type* mesh_n_element,
7473
const geo_element& K,
7474
std::vector<size_type>& i_dof) const;
7476
virtual bool is_continuous() const;
7477
virtual bool is_discontinuous() const { return !is_continuous(); }
7481
void dump(std::ostream& out = std::cerr) const;
7485
File: rheolef.info, Node: smart_pointer internal, Up: Internals
7487
8.17 `occurence', `smart_pointer' - memory management
7488
=====================================================
7490
(Source file: `util/lib/smart_pointer.h')
7495
Here is a convenient way to implement a true copy semantc, by using
7496
shallow copies and reference counting, in order to minimise memory
7497
copies. This concept is generally related to the "smart pointer"
7498
method for managing memory.
7500
The true semantic copy is defined as follows: if an object `A' is
7501
assigned to `B', such as `A = B', every further modification on `A' or
7502
`B' does not modify the other.
7508
class smart_pointer {
7513
smart_pointer (T* p = 0);
7514
smart_pointer (const smart_pointer&);
7516
smart_pointer& operator= (const smart_pointer&);
7520
const T* pointer () const;
7521
const T& data () const;
7522
const T* operator-> () const;
7523
const T& operator* () const;
7547
File: rheolef.info, Node: heap_allocator internal, Up: Internals
7549
8.18 heap_allocator - heap-based allocator
7550
==========================================
7552
(Source file: `util/lib/heap_allocator.h')
7557
Heap allocators are generally used when there is a lot of allocation
7558
and deallocation of small objects. For instance, this is often the
7559
case when dealing with `std::list' and `std::map'.
7561
Heap-based allocator is conform to the STL specification of allocators.
7562
It does not "free" the memory until the heap is destroyed.
7564
This allocator handles an a priori unlimited area of memory: a sequence
7565
of growing chunks are allocated. For a limited memory handler in
7566
the same spirit, see "stack_allocator"(9).
7571
typedef map <size_t, double, less<size_t>, heap_allocator<pair<size_t,double> > > map_type;
7573
a.insert (make_pair (0, 3.14));
7574
a.insert (make_pair (1, 1.17));
7575
for (map_type::iterator iter = a.begin(), last = a.end(); iter != last; iter++) {
7576
cout << (*iter).first << " " << (*iter).second << endl;
7582
template <typename T>
7583
class heap_allocator {
7585
struct handler_type; // forward declaration:
7590
typedef size_t size_type;
7591
typedef ptrdiff_t difference_type;
7593
typedef const T* const_pointer;
7594
typedef T& reference;
7595
typedef const T& const_reference;
7596
typedef T value_type;
7600
heap_allocator() throw()
7601
: handler (new handler_type)
7604
heap_allocator (const heap_allocator& ha) throw()
7605
: handler (ha.handler)
7607
++handler->reference_count;
7609
template <typename U>
7610
heap_allocator (const heap_allocator<U>& ha) throw()
7611
: handler ((typename heap_allocator<T>::handler_type*)(ha.handler))
7613
++handler->reference_count;
7615
~heap_allocator() throw()
7617
check_macro (handler != NULL, "unexpected null mem_info");
7618
if (--handler->reference_count == 0) delete handler;
7620
// Rebind to allocators of other types
7621
template <typename U>
7623
typedef heap_allocator<U> other;
7628
heap_allocator& operator= (const heap_allocator& ha)
7630
handler = ha.handler;
7631
++handler->reference_count;
7635
// utility functions:
7637
pointer address (reference r) const { return &r; }
7638
const_pointer address (const_reference c) const { return &c; }
7639
size_type max_size() const { return std::numeric_limits<size_t>::max() / sizeof(T); }
7641
// in-place construction/destruction
7643
void construct (pointer p, const_reference c)
7645
// placement new operator:
7646
new( reinterpret_cast<void*>(p) ) T(c);
7648
void destroy (pointer p)
7650
// call destructor directly:
7654
// allocate raw memory
7656
pointer allocate (size_type n, const void* = NULL)
7658
return pointer (handler->raw_allocate (n*sizeof(T)));
7660
void deallocate (pointer p, size_type n)
7662
// No need to free heap memory
7664
const handler_type* get_handler() const {
7671
handler_type* handler;
7672
template <typename U> friend class heap_allocator;
7676
File: rheolef.info, Node: stack_allocator internal, Up: Internals
7678
8.19 stack_allocator - stack-based allocator
7679
============================================
7681
(Source file: `util/lib/stack_allocator.h')
7686
Stack-based allocator, conform to the STL specification of allocators.
7687
Designed to use stack-based data passed as a parameter to the
7688
allocator constructor. Does not "free" the memory. Assumes that
7689
if the allocator is copied, stack memory is cleared and new allocations
7690
begin at the bottom of the stack again.
7692
Also works with any memory buffer, including heap memory. If the caller
7693
passes in heap memory, the caller is responsible for freeing the
7696
This allocator handles a limited area of memory: if this limit is
7697
reached, a "std::bad_alloc" exception is emmited. For a
7698
non-limited memory handler in the same spirit, see "heap_allocator"(9).
7703
const size_t stack_size = 1024;
7704
vector<unsigned char> stack (stack_size);
7705
stack_allocator<double> stack_alloc (stack.begin().operator->(), stack.size());
7706
typedef map <size_t, double, less<size_t>, stack_allocator<pair<size_t,double> > > map_type;
7707
map_type a (less<size_t>(), stack_alloc);
7708
a.insert (make_pair (0, 3.14));
7709
a.insert (make_pair (1, 1.17));
7710
for (map_type::iterator iter = a.begin(), last = a.end(); iter != last; iter++) {
7711
cout << (*iter).first << " " << (*iter).second << endl;
7717
template <typename T>
7718
class stack_allocator {
7720
struct handler_type; // forward declaration:
7725
typedef size_t size_type;
7726
typedef ptrdiff_t difference_type;
7728
typedef const T* const_pointer;
7729
typedef T& reference;
7730
typedef const T& const_reference;
7731
typedef T value_type;
7735
stack_allocator() throw()
7736
: handler (new handler_type)
7739
stack_allocator (unsigned char* stack, size_t stack_size) throw()
7740
: handler (new handler_type (stack, stack_size))
7742
warning_macro ("stack_allocator cstor");
7744
stack_allocator (const stack_allocator& sa) throw()
7745
: handler (sa.handler)
7747
++handler->reference_count;
7749
template <typename U>
7750
stack_allocator (const stack_allocator<U>& sa) throw()
7751
: handler ((typename stack_allocator<T>::handler_type*)(sa.handler))
7753
++handler->reference_count;
7755
~stack_allocator() throw()
7757
warning_macro ("stack_allocator dstor");
7758
check_macro (handler != NULL, "unexpected null mem_info");
7759
if (--handler->reference_count == 0) delete handler;
7761
// Rebind to allocators of other types
7762
template <typename U>
7764
typedef stack_allocator<U> other;
7769
stack_allocator& operator= (const stack_allocator& sa)
7771
handler = sa.handler;
7772
++handler->reference_count;
7776
// utility functions:
7778
pointer address (reference r) const { return &r; }
7779
const_pointer address (const_reference c) const { return &c; }
7780
size_type max_size() const { return std::numeric_limits<size_t>::max() / sizeof(T); }
7782
// in-place construction/destruction
7784
void construct (pointer p, const_reference c)
7786
// placement new operator:
7787
new( reinterpret_cast<void*>(p) ) T(c);
7789
void destroy (pointer p)
7791
// call destructor directly:
7795
// allocate raw memory
7797
pointer allocate (size_type n, const void* = NULL)
7799
warning_macro ("allocate "<<n<<" type " << typename_macro(T));
7800
check_macro (handler->stack != NULL, "unexpected null stack");
7801
void* p = handler->stack + handler->allocated_size;
7802
handler->allocated_size += n*sizeof(T);
7804
if (handler->allocated_size + 1 > handler->max_size) {
7805
warning_macro ("stack is full: throwing...");
7806
throw std::bad_alloc();
7810
void deallocate (pointer p, size_type n)
7812
warning_macro ("deallocate "<<n<<" type "<<typename_macro(T));
7813
// No need to free stack memory
7815
const handler_type* get_handler() const {
7822
struct handler_type {
7823
unsigned char* stack;
7824
size_t allocated_size;
7826
size_t reference_count;
7834
warning_macro ("stack_allocator::mem_info cstor NULL");
7836
handler_type (unsigned char* stack1, size_t size1)
7842
warning_macro ("stack_allocator::mem_info cstori: size="<<max_size);
7846
warning_macro ("stack_allocator::mem_info dstor: size="<<max_size);
7849
handler_type* handler;
7850
template <typename U> friend class stack_allocator;
7853
template <typename T1>
7854
bool operator==( const stack_allocator<T1>& lhs, const stack_allocator<T1>& rhs) throw()
7856
return lhs.get_handler() == rhs.get_handler();
7858
template <typename T1>
7859
bool operator!=( const stack_allocator<T1>& lhs, const stack_allocator<T1>& rhs) throw()
7861
return lhs.get_handler() != rhs.get_handler();
7865
File: rheolef.info, Node: pretty_name internal, Up: Internals
7867
8.20 `typename_macro', `pretty_typename_macro' - type demangler and pretty printer
7868
==================================================================================
7870
(Source file: `util/lib/pretty_name.h')
7875
These preprocessor macro-definitions are usefull when dealing with
7876
complex types as generated by imbricted template technics: they print
7877
in clear a complex type at run-time. `typeid_name_macro' obtains a
7878
human readable type in a `std::tring' form by calling the system
7879
`typeid' function and then a demangler. When this type is very long,
7880
`pretty_name_macro' prints also it in a multi-line form with a pretty
7886
typedef map <size_t, double, less<size_t>, heap_allocator<pair<size_t,double> > > map_type;
7887
cout << typeid_name_macro (map_type);
7892
extern std::string typeid_name (const char* name, bool do_indent);
7894
/// @brief get string from a type, with an optional pretty-printing for complex types
7895
#define typename_macro(T) typeid_name(typeid(T).name(), false)
7896
#define pretty_typename_macro(T) typeid_name(typeid(T).name(), true)
7898
/// @brief get string type from a variable or expression
7899
template <class T> std::string typename_of (T x) { return typename_macro(T); }
7900
template <class T> std::string pretty_typename_of (T x) { return pretty_typename_macro(T); }
7903
File: rheolef.info, Node: acinclude internal, Up: Internals
7905
8.21 `acinclude' - autoconf macros
7906
==================================
7908
(Source file: `config/acinclude.m4')
7913
These macros test for particular system featutres that rheolef
7914
uses. These tests print the messages telling the user which feature
7915
they are looking for and what they find. They cache their results
7916
for future `configure' runs. Some of these macros "set"
7917
some shell variable, "defines" some output variables for the
7918
`config.h' header, or performs Makefile macros "subsitutions".
7919
See `autoconf' documentation for how to use such variables.
7924
Follows a list of particular check required for a successfull
7929
Check to see if GiNaC library exists. If so, set the shell
7930
variable `rheo_have_ginac' to "yes", defines HAVE_GINAC and
7931
substitues INCLUDES_GINAC and LADD_GINAC for adding in CFLAGS and
7932
LDFLAGS, respectively, If not, set the shell variable
7933
rheo_have_ginac to "no".
7937
Check to see if library `-lcln' exists. If so, set the shell
7938
variable `rheo_have_cln' to "yes", defines HAVE_CLN and
7939
substitues INCLUDES_CLN and LADD_CLN for adding in CFLAGS and
7940
LDFLAGS, respectively, If not, set the shell variable no "no".
7941
Includes and libraries path are searched from a given shell
7942
variable `rheo_dir_cln'. This shell variable could be set for
7943
instance by an appropriate `--with-cln'=VALUE_DIR_CLN option.
7944
The default value is `/usr/local/math'.
7946
`RHEO_CHECK_SPOOLES_2_0'
7948
Check to see if spooles library has old version 2.0 since
7949
`FrontMtx_factorInpMtx' profile has changed in version 2.2. If
7950
so, defines HAVE_SPOOLES_2_0. This macro is called by
7955
Check to see if taucs library and headers exists. If so, set the
7956
shell variable "rheo_have_taucs" to "yes", defines HAVE_TAUCS and
7957
substitues INCLUDES_TAUCS and LADD_TAUCS for adding in CXXFLAGS
7958
and LDFLAGS, respectively, If not, set the shell variable to "no".
7959
Includes and libraries options are given shell variable
7960
$rheo_ldadd_taucs and $rheo_incdir_taucs. These shell variables
7961
could be set for instance by appropriates
7962
"-with-taucs-ldadd="'rheo_ldadd_taucs' and
7963
"-with-taucs-includes="'rheo_incdir_taucs' options.
7967
Check to see if boost headers exists. If so, set the shell
7968
variable "rheo_have_boost" to "yes", defines HAVE_BOOST and
7969
substitues INCLUDES_BOOST for adding in CXXFLAGS, and LDADD_BOOST
7970
for adding in LIBS. If not, set the shell variable to "no".
7971
Includes options are given in the shell variables
7972
$rheo_incdir_boost and $rheo_libdir_boost. These shell variables
7973
could be set for instance by appropriates
7974
"-with-boost-includes="'rheo_incdir_boost' and
7975
"-with-boost-libdir="'rheo_libdir_boost' options.
7979
Check to see if zlib library and headers exists. If so, set the
7980
shell variable "rheo_have_zlib" to "yes", defines HAVE_ZLIB and
7981
substitues INCLUDES_ZLIB and LADD_ZLIB for adding in CXXFLAGS and
7982
LDFLAGS, respectively, If not, set the shell variable to "no".
7983
Includes and libraries path are searched from given shell variable
7984
$rheo_dir_zlib/lib and $rheo_incdir_zlib. Default value for
7985
$rheo_incdir_zlib is $rheo_dir_zlib/include. These shell variables
7986
could be set for instance by appropriates "-with-zlib="'dir_zlib'
7987
and "-with-zlib-includes="'incdir_zlib' options.
7989
`RHEO_CHECK_SPOOLES'
7991
Check to see if spooles library and headers exists. If so, set the
7992
shell variable "rheo_have_spooles" to "yes", defines HAVE_SPOOLES
7993
and substitues INCLUDES_SPOOLES and LADD_SPOOLES for adding in
7994
CXXFLAGS and LDFLAGS, respectively, If not, set the shell variable
7995
to "no". Includes and libraries path are searched from given
7996
shell variable "rheo_libdir_spooles" and "rheo_incdir_spooles".
7997
These shell variables could be set for instance by appropriates
7998
"-with-spooles="'libdir_spooles' and
7999
"-with-spooles-includes="'incdir_spooles' options.
8001
`RHEO_CHECK_UMFPACK'
8003
Check to see if umfpack library and headers exists. If so, set the
8004
shell variable "rheo_have_umfpack" to "yes", defines HAVE_UMFPACK
8005
and substitues INCLUDES_UMFPACK and LADD_UMFPACK for adding in
8006
CXXFLAGS and LDFLAGS, respectively, If not, set the shell variable
8007
to "no". Includes and libraries path are searched from given
8008
shell variable "rheo_libdir_umfpack" and "rheo_incdir_umfpack".
8009
These shell variables could be set for instance by appropriates
8010
"-with-umfpack="'libdir_umfpack' and
8011
"-with-umfpack-includes="'incdir_umfpack' options.
8013
`RHEO_CHECK_MALLOC_DBG'
8015
Check to see if malloc debug library -lmalloc_dbg and corresponding
8016
header <malloc_dbg.h> exists. If so, set the shell variable
8017
`rheo_have_malloc_dbg' to "yes", defines HAVE_MALLOC_DBG, add
8018
`-I'DIR_MALLOC_DBG`/include' to CFLAGS, add
8019
DIR_MALLOC_DBG`/lib/libmalloc_dbg.a' to LDFLAGS. Here,
8020
DIR_MALLOC_DBG is the directory such that DIR_MALLOC_DBG`/bin'
8021
appears in PATH and the command DIR_MALLOC_DBG`/bin/malloc_dbg'
8022
exists. If not, set the variable to "no". Set also
8023
LIBS_MALLOC_DBG to these flags.
8025
`RHEO_CHECK_DMALLOC'
8027
Check whether the dmalloc package exists and set the corresponding
8028
shell value "rheo_have_dmalloc" and HAVE_DMALLOC (in Makefile.am
8029
and config.h) accordingly, create LDADD_DMALLOC and LDADD_DMALLOCXX
8030
Makefile.am variables.
8032
`RHEO_CHECK_NAMESPACE'
8034
Check whether the namespace feature is supported by the C++
8035
compiler. value. So, try to compile the following code:
8036
namespace computers {
8037
struct keyboard { int getkey() const { return 0; } };
8040
struct keyboard { void playNote(int note); };
8043
void keyboard::playNote(int note) { }
8045
using namespace computers;
8053
If it compile, set the corresponding shell variable
8054
"rheo_have_namespace" to "yes" and defines HAVE_NAMESPACE. If
8055
not, set the variable no "no".
8057
`RHEO_CHECK_STD_NAMESPACE'
8059
Some compilers (e.g. GNU C++ 2.7.2) does not support the full
8060
namespace feature. Nevertheless, they support the "std:" namespace
8061
for the C++ library. is supported by the C++ compiler. The
8062
following code is submitted to the compiler:
8064
extern "C" void exit(int);
8066
std::vector<int> x(3);
8069
If it compile, set the corresponding shell variable
8070
"rheo_have_std_namespace" to "yes" and defines HAVE_STD_NAMESPACE.
8071
If not, set the variable no "no".
8073
`RHEO_PROG_GNU_MAKE'
8075
Find command make and check whether make is GNU make. If so, set
8076
the corresponding shell variable "rheo_prog_gnu_make" to "yes"
8077
and substitues no_print_directory_option to "-no-print-directory".
8078
If not, set the shell variable no "no".
8080
`RHEO_CHECK_ISTREAM_RDBUF'
8084
Check to see if "iostream::rdbuf(void*)" function exists, that set
8085
the "ios" buffer of a stream. Despite this function is standard,
8086
numerous compilers does not furnish it. a common implementation is
8087
to set directly the buffer variable. For instance, the CRAY C++
8088
compiler implements this variable as "ios::bp". These two
8089
functions set the shell variables "rheo_have_istream_rdbuf" and
8090
"rheo_have_ios_bp" and define HAVE_ISTREAM_RDBUF and HAVE_IOS_BP
8093
`RHEO_CHECK_IOS_SETSTATE'
8095
Check to see if "ios::setstate(long)" function exists, that set the
8096
"ios" state variable of a stream. Despite this function is standard,
8097
numerous compilers does not furnish it. a common implementation
8098
is to set directly the buffer variable. For instance, the CRAY C++
8099
compiler does not implements it. This function set the shell
8100
variables "rheo_have_ios_setstate" and define HAVE_IOS_SETSTATE.
8102
`RHEO_CHECK_FILEBUF_INT'
8104
`RHEO_CHECK_FILEBUF_FILE'
8106
`RHEO_CHECK_FILEBUF_FILE_MODE'
8108
Check wheter "filebuf::filebuf(int fileno)",
8109
"filebuf::filebuf(FILE* fd)" exist, or "filebuf::filebuf(FILE* fd,
8110
ios::openmode)" exist, respectively. If so, set the
8111
corresponding shell variable "rheo_have_filebuf_int" (resp.
8112
"rheo_have_filebuf_file") to "yes" and defines HAVE_FILEBUF_INT,
8113
(resp. HAVE_FILEBUF_FILE). If not, set the variable no "no".
8114
Notes that there is no standardisation of this function in the
8115
"c++" library. Nevertheless, this fonctionality is usefull to
8116
open a pipe stream class, as "pstream(3)".
8118
`RHEO_CHECK_GETTIMEOFDAY'
8120
Check whether the "gettimeofday(timeval*, timezone*)" function
8121
exists and set the corresponding shell value
8122
"rheo_have_gettimeofday" and define HAVE_GETTIMEOFDAY accordingly.
8124
`RHEO_CHECK_WIERDGETTIMEOFDAY'
8126
This is for Solaris, where they decided to change the CALLING
8127
SEQUENCE OF gettimeofday! Check whether the
8128
"gettimeofday(timeval*)" function exists and set the corresponding
8129
shell value "rheo_have_wierdgettimeofday" and define
8130
HAVE_WIERDGETTIMEOFDAY accordingly.
8132
`RHEO_CHECK_BSDGETTIMEOFDAY'
8134
For BSD systems, check whether the "BSDgettimeofday(timeval*,
8135
timezone*)" function exists and set the corresponding shell
8136
value "rheo_have_bsdgettimeofday" and define HAVE_BSDGETTIMEOFDAY
8139
`RHEO_CHECK_AMICCLK'
8141
Check whether the clock "amicclk()" function exists and set the
8142
corresponding shell value "rheo_have_amicclk" and define
8143
HAVE_AMICCLK accordingly.
8145
`RHEO_CHECK_TEMPLATE_FULL_SPECIALIZATION'
8147
Check whether the template specialization syntax "template<>" is
8148
supported by the compiler value. So, try to compile, run and
8149
check the return value for the following code:
8150
template<class T> struct toto {
8151
int tutu() const { return 1; }
8153
template<> struct toto<float> {
8154
int tutu() const { return 0; }
8160
If so, set the corresponding shell variable
8161
"rheo_have_template_full_specialization" to "yes" and defines
8162
HAVE_TEMPLATE_FULL_SPECIALIZATION. If not, set the variable no
8165
`RHEO_CHECK_ISNAN_DOUBLE'
8167
`RHEO_CHECK_ISINF_DOUBLE'
8169
`RHEO_CHECK_FINITE_DOUBLE'
8171
`RHEO_CHECK_INFINITY'
8173
`RHEO_CHECK_ABS_DOUBLE'
8175
`RHEO_CHECK_SQR_DOUBLE'
8177
Check whether the funtions
8180
bool finite(double);
8184
are supported by the compiler, respectively. If so, set the
8185
corresponding shell variable "rheo_have_xxx" to "yes" and
8186
defines HAVE_XXX. If not, set the variable no "no".
8190
Check to see if the "flex" command and the corresponding header
8191
"FlexLexer.h" are available. If so, set the shell variable
8192
"rheo_have_flex" to "yes" and substitues FLEX to "flex" and
8193
FLEXLEXER_H to the full path for FlexLexer.h If not, set the shell
8198
Check wheter we are using KAI C++ compiler. If so, set the shell
8199
variable "ac_cv_prog_kcc" to "yes". If not, set the shell
8200
variable no "no". The shell variable is also exported for
8201
sub-shells, such as ltconfig from libtool.
8205
Check wheter we are using CRAY C++ compiler. If so, set the shell
8206
variable "ac_cv_prog_cray_cc" to "yes" and defines HAVE_CRAY_CXX.
8207
If not, set the shell variable no "no". The shell variable is
8208
also exported for sub-shells.
8212
Check wheter we are using DEC C++ compiler. If so, set the shell
8213
variable "ac_cv_prog_dec_cc" to "yes". If not, set the shell
8214
variable no "no". The shell variable is also exported for
8215
sub-shells, such as ltconfig from libtool.
8217
`RHEO_RECOGNIZE_CXX'
8219
Check wheter we are able to recognize the C++ compiler. Tested
8221
The KAI C++ compiler that defines: __KCC (KCC-3.2b)
8222
The GNU C++ compiler that defines: __GNUC__ (egcs-1.1.1)
8223
The CRAY C++ compiler that defines: cray
8224
The DEC C++ compiler that defines: __DECCXX
8226
If so, substitue RECOGNIZED_CXX to a specific compiler's rule file,
8227
e.g, "${top_srcdir}/config/gnu_cxx.mk", for a subsequent Makefile
8228
include. If not, substitue to "/dev/null". Substitutes also
8229
EXTRA_LDFLAGS. Raw cc is the C compiler associated to the C++ one.
8230
By this way C and C++ files are handled with a .c suffix. Special C
8231
files that requiere the cc compiler, such as "alloca.c" use some
8232
specific makefile rule.
8235
AC_PROG_CC(gcc cc cl)
8236
AC_PROG_CXX(c++ g++ cxx KCC CC CC cc++ xlC aCC)
8241
Set some optimization flags associated to the recognized C++ compiler
8244
`RHEO_CHECK_LATEX_HYPEREF'
8246
Check whether the hyperref LaTeX package exists and set the
8247
corresponding shell value "rheo_have_latex_hyperref" and
8248
HAVE_LATEX_HYPEREF (for Makefiles) accordingly.
8252
Check for the "mpirun" command, the corresponding header "mpi.h"
8253
and library "-lmpi" are available. If so, set the shell variable
8254
"rheo_have_mpi" to "yes", defines HAVE_MPI and substitues MPIRUN to
8255
"mpirun" and RUN to "mpirun -np 2", INCLUDES_MPI and LDADD_MPI.
8256
If not, set the shell variable no "no".
8258
`RHEO_CHECK_PARMETIS'
8260
Check for the "parmetis" and "metis" libraries. Defines
8261
HAVE_PARMETIS and substitues INCLUDES_MPI and LDADD_MPI.
8262
Requires the MPI library.
8266
Check for the "scotch" distributed mesh partitionner libraries.
8267
Defines HAVE_SCOTCH and substitues INCLUDES_SCOTCH and LDADD_SCOTCH.
8268
Requires the MPI library.
8272
Check for the "blas" basic linear algebra subroutines library.
8273
Defines HAVE_BLAS and substitues LDADD_BLAS and INCLUDES_BLAS.
8277
Check for the "pastix" distributed direct solver libraries.
8278
Defines HAVE_PASTIX and substitues INCLUDES_PASTIX and LDADD_PASTIX.
8279
Requires the MPI and SCOTCH libraries.
8282
File: rheolef.info, Node: FAQ for developers, Up: Top
8284
9 FAQ for developers
8285
********************
8287
This list of Frequently Asked Questions intended for Rheolef developers
8288
and maintainers. I'm looking for new questions (_with_ answers, better
8289
answers, or both. Please, send suggestions to
8290
<Pierre.Saramito@imag.fr>.
8292
9.1 How to regenerate the `configure' script
8293
============================================
8295
The configure script and makefiles are automatically produced from file
8296
`configure.ac' and `Makefile.am' by using the autoconf and automake
8301
9.1.1 In which order does the things build ?
8302
--------------------------------------------
8304
Let us look at with details the configure files flow:
8306
[acinclude.m4] -----> aclocal* -----------> [aclocal.m4]
8309
----------------+---> autoconf* ----------> configure
8312
[Makefile.am] ------> automake* ----------> Makefile.in
8314
[config.h.in] -+ +-> config.h
8316
Makefile.in ---+----> configure* -------+-> Makefile
8318
[config.mk.in] + +-> config.mk
8320
9.1.2 What means these commands ?
8321
---------------------------------
8323
Let us review the list of commands:
8326
take `acinclude.m4' and build `aclocal.m4'. The
8327
arguments specifies that these files are located in the
8328
`config/' directory.
8331
translate every `Makefile.am' and then build a `Makefile.in'.
8334
translate `configure.ac' in `configure' which is an
8335
executable shell script. The arguments specifies that
8336
`aclocal.m4' is located in the `config/' dirctory.
8337
All this files are machine independent.
8340
the automatically generated script, scan the machine and
8341
translate `config.h.in', `config.mk.in' and every
8342
`Makefile.in' into `config.h', `config.mk' and every
8343
`Makefile', respectively. At this step, all produced
8344
files are machine dependent.
8347
9.2 How to save my version ?
8348
============================
8350
First, check that our distribution is valid
8353
First, check that your modifications are not in conflict with others.
8354
Go to the top source directory and enter:
8358
9.2.1 Easy: no conflicts with another developer
8359
-----------------------------------------------
8361
A listing of labeled files appears:
8363
Modified skit/lib/blas1_tst.c
8364
Update skit/lib/blas2_tst.c
8366
Its mean that you have modified `blas1_tst.c'. Another concurrent
8367
developer has modified `blas2_tst.c', and your local file version is
8368
not up-to-date. There is no conflict, labeled by a `*Merge*' label.
8370
First, update your local version:
8374
Before to store your version of the Rheolef distribution, check the
8375
consistency by running rnon-regression tests:
8379
When all tests are ok:
8383
and enter a change log comment terminated by the `ctrl-D' key.
8385
Check now that your version status is the most up-to-date Rheolef
8390
9.2.2 I have conflicts with another developer
8391
---------------------------------------------
8393
The listing issued by `make status' looks like:
8395
Modified skit/lib/blas1_tst.c
8396
Update skit/lib/blas2_tst.c
8397
*Merge* skit/lib/blas3_tst.c
8399
Its mean that you and another developer have modified at least one
8400
common line in `blas3_tst.c'. Moreover, the developer has already
8401
saved his version in a previous Rheolef distribution. You have now to
8402
merge these modifications. Thus, enter:
8404
mv blas3_tst.c blas3_tst.c.new
8405
cvs checkout blas3_tst.c
8406
sdiff blas3_tst.c blas3_tst.c.new | more
8407
and then edit `blas3_tst.c' to integrate the modifications of
8408
`blas3_tst.c.new', generated by another developer. When it is done, go
8409
to the top directory and enter,
8412
Modified skit/lib/blas1_tst.c
8413
Update skit/lib/blas2_tst.c
8414
Modified skit/lib/blas3_tst.c
8415
The situation becomes mature:
8417
It will update the `blas2_tst.c'. Finally,
8419
that will save modified `blas1_tst.c' and `blas3_tst.c'.
8421
9.2.3 I have deleted a source file...
8422
-------------------------------------
8430
How to restaure the file, now ?
8434
cvs checkout Makefile.am
8436
You restaure the last available version of the missing file in the most
8437
up-to-date Rheolef distribution.
8440
File: rheolef.info, Node: Copying, Up: Top
8442
Annexe A GNU General Public License
8443
***********************************
8445
Version 2, June 1991
8447
Copyright (C) 1989, 1991 Free Software Foundation, Inc.
8448
675 Mass Ave, Cambridge, MA 02139, USA
8450
Everyone is permitted to copy and distribute verbatim copies
8451
of this license document, but changing it is not allowed.
8456
The licenses for most software are designed to take away your freedom
8457
to share and change it. By contrast, the GNU General Public License is
8458
intended to guarantee your freedom to share and change free
8459
software--to make sure the software is free for all its users. This
8460
General Public License applies to most of the Free Software
8461
Foundation's software and to any other program whose authors commit to
8462
using it. (Some other Free Software Foundation software is covered by
8463
the GNU Library General Public License instead.) You can apply it to
8466
When we speak of free software, we are referring to freedom, not price.
8467
Our General Public Licenses are designed to make sure that you have the
8468
freedom to distribute copies of free software (and charge for this
8469
service if you wish), that you receive source code or can get it if you
8470
want it, that you can change the software or use pieces of it in new
8471
free programs; and that you know you can do these things.
8473
To protect your rights, we need to make restrictions that forbid anyone
8474
to deny you these rights or to ask you to surrender the rights. These
8475
restrictions translate to certain responsibilities for you if you
8476
distribute copies of the software, or if you modify it.
8478
For example, if you distribute copies of such a program, whether gratis
8479
or for a fee, you must give the recipients all the rights that you
8480
have. You must make sure that they, too, receive or can get the source
8481
code. And you must show them these terms so they know their rights.
8483
We protect your rights with two steps: (1) copyright the software, and
8484
(2) offer you this license which gives you legal permission to copy,
8485
distribute and/or modify the software.
8487
Also, for each author's protection and ours, we want to make certain
8488
that everyone understands that there is no warranty for this free
8489
software. If the software is modified by someone else and passed on, we
8490
want its recipients to know that what they have is not the original, so
8491
that any problems introduced by others will not reflect on the original
8492
authors' reputations.
8494
Finally, any free program is threatened constantly by software patents.
8495
We wish to avoid the danger that redistributors of a free program will
8496
individually obtain patent licenses, in effect making the program
8497
proprietary. To prevent this, we have made it clear that any patent
8498
must be licensed for everyone's free use or not licensed at all.
8500
The precise terms and conditions for copying, distribution and
8501
modification follow.
8503
TERMS AND CONDITIONS FOR COPYING, DISTRIBUTION AND MODIFICATION
8504
0. This License applies to any program or other work which contains a
8505
notice placed by the copyright holder saying it may be distributed
8506
under the terms of this General Public License. The "Program",
8507
below, refers to any such program or work, and a "work based on
8508
the Program" means either the Program or any derivative work under
8509
copyright law: that is to say, a work containing the Program or a
8510
portion of it, either verbatim or with modifications and/or
8511
translated into another language. (Hereinafter, translation is
8512
included without limitation in the term "modification".) Each
8513
licensee is addressed as "you".
8515
Activities other than copying, distribution and modification are
8516
not covered by this License; they are outside its scope. The act
8517
of running the Program is not restricted, and the output from the
8518
Program is covered only if its contents constitute a work based on
8519
the Program (independent of having been made by running the
8520
Program). Whether that is true depends on what the Program does.
8522
1. You may copy and distribute verbatim copies of the Program's
8523
source code as you receive it, in any medium, provided that you
8524
conspicuously and appropriately publish on each copy an appropriate
8525
copyright notice and disclaimer of warranty; keep intact all the
8526
notices that refer to this License and to the absence of any
8527
warranty; and give any other recipients of the Program a copy of
8528
this License along with the Program.
8530
You may charge a fee for the physical act of transferring a copy,
8531
and you may at your option offer warranty protection in exchange
8534
2. You may modify your copy or copies of the Program or any portion
8535
of it, thus forming a work based on the Program, and copy and
8536
distribute such modifications or work under the terms of Section 1
8537
above, provided that you also meet all of these conditions:
8539
a. You must cause the modified files to carry prominent notices
8540
stating that you changed the files and the date of any change.
8542
b. You must cause any work that you distribute or publish, that
8543
in whole or in part contains or is derived from the Program
8544
or any part thereof, to be licensed as a whole at no charge
8545
to all third parties under the terms of this License.
8547
c. If the modified program normally reads commands interactively
8548
when run, you must cause it, when started running for such
8549
interactive use in the most ordinary way, to print or display
8550
an announcement including an appropriate copyright notice and
8551
a notice that there is no warranty (or else, saying that you
8552
provide a warranty) and that users may redistribute the
8553
program under these conditions, and telling the user how to
8554
view a copy of this License. (Exception: if the Program
8555
itself is interactive but does not normally print such an
8556
announcement, your work based on the Program is not required
8557
to print an announcement.)
8559
These requirements apply to the modified work as a whole. If
8560
identifiable sections of that work are not derived from the
8561
Program, and can be reasonably considered independent and separate
8562
works in themselves, then this License, and its terms, do not
8563
apply to those sections when you distribute them as separate
8564
works. But when you distribute the same sections as part of a
8565
whole which is a work based on the Program, the distribution of
8566
the whole must be on the terms of this License, whose permissions
8567
for other licensees extend to the entire whole, and thus to each
8568
and every part regardless of who wrote it.
8570
Thus, it is not the intent of this section to claim rights or
8571
contest your rights to work written entirely by you; rather, the
8572
intent is to exercise the right to control the distribution of
8573
derivative or collective works based on the Program.
8575
In addition, mere aggregation of another work not based on the
8576
Program with the Program (or with a work based on the Program) on
8577
a volume of a storage or distribution medium does not bring the
8578
other work under the scope of this License.
8580
3. You may copy and distribute the Program (or a work based on it,
8581
under Section 2) in object code or executable form under the terms
8582
of Sections 1 and 2 above provided that you also do one of the
8585
a. Accompany it with the complete corresponding machine-readable
8586
source code, which must be distributed under the terms of
8587
Sections 1 and 2 above on a medium customarily used for
8588
software interchange; or,
8590
b. Accompany it with a written offer, valid for at least three
8591
years, to give any third party, for a charge no more than your
8592
cost of physically performing source distribution, a complete
8593
machine-readable copy of the corresponding source code, to be
8594
distributed under the terms of Sections 1 and 2 above on a
8595
medium customarily used for software interchange; or,
8597
c. Accompany it with the information you received as to the offer
8598
to distribute corresponding source code. (This alternative is
8599
allowed only for noncommercial distribution and only if you
8600
received the program in object code or executable form with
8601
such an offer, in accord with Subsection b above.)
8603
The source code for a work means the preferred form of the work for
8604
making modifications to it. For an executable work, complete
8605
source code means all the source code for all modules it contains,
8606
plus any associated interface definition files, plus the scripts
8607
used to control compilation and installation of the executable.
8608
However, as a special exception, the source code distributed need
8609
not include anything that is normally distributed (in either
8610
source or binary form) with the major components (compiler,
8611
kernel, and so on) of the operating system on which the executable
8612
runs, unless that component itself accompanies the executable.
8614
If distribution of executable or object code is made by offering
8615
access to copy from a designated place, then offering equivalent
8616
access to copy the source code from the same place counts as
8617
distribution of the source code, even though third parties are not
8618
compelled to copy the source along with the object code.
8620
4. You may not copy, modify, sublicense, or distribute the Program
8621
except as expressly provided under this License. Any attempt
8622
otherwise to copy, modify, sublicense or distribute the Program is
8623
void, and will automatically terminate your rights under this
8624
License. However, parties who have received copies, or rights,
8625
from you under this License will not have their licenses
8626
terminated so long as such parties remain in full compliance.
8628
5. You are not required to accept this License, since you have not
8629
signed it. However, nothing else grants you permission to modify
8630
or distribute the Program or its derivative works. These actions
8631
are prohibited by law if you do not accept this License.
8632
Therefore, by modifying or distributing the Program (or any work
8633
based on the Program), you indicate your acceptance of this
8634
License to do so, and all its terms and conditions for copying,
8635
distributing or modifying the Program or works based on it.
8637
6. Each time you redistribute the Program (or any work based on the
8638
Program), the recipient automatically receives a license from the
8639
original licensor to copy, distribute or modify the Program
8640
subject to these terms and conditions. You may not impose any
8641
further restrictions on the recipients' exercise of the rights
8642
granted herein. You are not responsible for enforcing compliance
8643
by third parties to this License.
8645
7. If, as a consequence of a court judgment or allegation of patent
8646
infringement or for any other reason (not limited to patent
8647
issues), conditions are imposed on you (whether by court order,
8648
agreement or otherwise) that contradict the conditions of this
8649
License, they do not excuse you from the conditions of this
8650
License. If you cannot distribute so as to satisfy simultaneously
8651
your obligations under this License and any other pertinent
8652
obligations, then as a consequence you may not distribute the
8653
Program at all. For example, if a patent license would not permit
8654
royalty-free redistribution of the Program by all those who
8655
receive copies directly or indirectly through you, then the only
8656
way you could satisfy both it and this License would be to refrain
8657
entirely from distribution of the Program.
8659
If any portion of this section is held invalid or unenforceable
8660
under any particular circumstance, the balance of the section is
8661
intended to apply and the section as a whole is intended to apply
8662
in other circumstances.
8664
It is not the purpose of this section to induce you to infringe any
8665
patents or other property right claims or to contest validity of
8666
any such claims; this section has the sole purpose of protecting
8667
the integrity of the free software distribution system, which is
8668
implemented by public license practices. Many people have made
8669
generous contributions to the wide range of software distributed
8670
through that system in reliance on consistent application of that
8671
system; it is up to the author/donor to decide if he or she is
8672
willing to distribute software through any other system and a
8673
licensee cannot impose that choice.
8675
This section is intended to make thoroughly clear what is believed
8676
to be a consequence of the rest of this License.
8678
8. If the distribution and/or use of the Program is restricted in
8679
certain countries either by patents or by copyrighted interfaces,
8680
the original copyright holder who places the Program under this
8681
License may add an explicit geographical distribution limitation
8682
excluding those countries, so that distribution is permitted only
8683
in or among countries not thus excluded. In such case, this
8684
License incorporates the limitation as if written in the body of
8687
9. The Free Software Foundation may publish revised and/or new
8688
versions of the General Public License from time to time. Such
8689
new versions will be similar in spirit to the present version, but
8690
may differ in detail to address new problems or concerns.
8692
Each version is given a distinguishing version number. If the
8693
Program specifies a version number of this License which applies
8694
to it and "any later version", you have the option of following
8695
the terms and conditions either of that version or of any later
8696
version published by the Free Software Foundation. If the Program
8697
does not specify a version number of this License, you may choose
8698
any version ever published by the Free Software Foundation.
8700
10. If you wish to incorporate parts of the Program into other free
8701
programs whose distribution conditions are different, write to the
8702
author to ask for permission. For software which is copyrighted
8703
by the Free Software Foundation, write to the Free Software
8704
Foundation; we sometimes make exceptions for this. Our decision
8705
will be guided by the two goals of preserving the free status of
8706
all derivatives of our free software and of promoting the sharing
8707
and reuse of software generally.
8710
11. BECAUSE THE PROGRAM IS LICENSED FREE OF CHARGE, THERE IS NO
8711
WARRANTY FOR THE PROGRAM, TO THE EXTENT PERMITTED BY APPLICABLE
8712
LAW. EXCEPT WHEN OTHERWISE STATED IN WRITING THE COPYRIGHT
8713
HOLDERS AND/OR OTHER PARTIES PROVIDE THE PROGRAM "AS IS" WITHOUT
8714
WARRANTY OF ANY KIND, EITHER EXPRESSED OR IMPLIED, INCLUDING, BUT
8715
NOT LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND
8716
FITNESS FOR A PARTICULAR PURPOSE. THE ENTIRE RISK AS TO THE
8717
QUALITY AND PERFORMANCE OF THE PROGRAM IS WITH YOU. SHOULD THE
8718
PROGRAM PROVE DEFECTIVE, YOU ASSUME THE COST OF ALL NECESSARY
8719
SERVICING, REPAIR OR CORRECTION.
8721
12. IN NO EVENT UNLESS REQUIRED BY APPLICABLE LAW OR AGREED TO IN
8722
WRITING WILL ANY COPYRIGHT HOLDER, OR ANY OTHER PARTY WHO MAY
8723
MODIFY AND/OR REDISTRIBUTE THE PROGRAM AS PERMITTED ABOVE, BE
8724
LIABLE TO YOU FOR DAMAGES, INCLUDING ANY GENERAL, SPECIAL,
8725
INCIDENTAL OR CONSEQUENTIAL DAMAGES ARISING OUT OF THE USE OR
8726
INABILITY TO USE THE PROGRAM (INCLUDING BUT NOT LIMITED TO LOSS OF
8727
DATA OR DATA BEING RENDERED INACCURATE OR LOSSES SUSTAINED BY YOU
8728
OR THIRD PARTIES OR A FAILURE OF THE PROGRAM TO OPERATE WITH ANY
8729
OTHER PROGRAMS), EVEN IF SUCH HOLDER OR OTHER PARTY HAS BEEN
8730
ADVISED OF THE POSSIBILITY OF SUCH DAMAGES.
8732
END OF TERMS AND CONDITIONS
8733
How to Apply These Terms to Your New Programs
8734
=============================================
8736
If you develop a new program, and you want it to be of the greatest
8737
possible use to the public, the best way to achieve this is to make it
8738
free software which everyone can redistribute and change under these
8741
To do so, attach the following notices to the program. It is safest to
8742
attach them to the start of each source file to most effectively convey
8743
the exclusion of warranty; and each file should have at least the
8744
"copyright" line and a pointer to where the full notice is found.
8746
ONE LINE TO GIVE THE PROGRAM'S NAME AND AN IDEA OF WHAT IT DOES.
8747
Copyright (C) 19YY NAME OF AUTHOR
8749
This program is free software; you can redistribute it and/or
8750
modify it under the terms of the GNU General Public License
8751
as published by the Free Software Foundation; either version 2
8752
of the License, or (at your option) any later version.
8754
This program is distributed in the hope that it will be useful,
8755
but WITHOUT ANY WARRANTY; without even the implied warranty of
8756
MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
8757
GNU General Public License for more details.
8759
You should have received a copy of the GNU General Public License
8760
along with this program; if not, write to the Free Software
8761
Foundation, Inc., 675 Mass Ave, Cambridge, MA 02139, USA.
8763
Also add information on how to contact you by electronic and paper mail.
8765
If the program is interactive, make it output a short notice like this
8766
when it starts in an interactive mode:
8768
Gnomovision version 69, Copyright (C) 19YY NAME OF AUTHOR
8769
Gnomovision comes with ABSOLUTELY NO WARRANTY; for details
8770
type `show w'. This is free software, and you are welcome
8771
to redistribute it under certain conditions; type `show c'
8774
The hypothetical commands `show w' and `show c' should show the
8775
appropriate parts of the General Public License. Of course, the
8776
commands you use may be called something other than `show w' and `show
8777
c'; they could even be mouse-clicks or menu items--whatever suits your
8780
You should also get your employer (if you work as a programmer) or your
8781
school, if any, to sign a "copyright disclaimer" for the program, if
8782
necessary. Here is a sample; alter the names:
8784
Yoyodyne, Inc., hereby disclaims all copyright
8785
interest in the program `Gnomovision'
8786
(which makes passes at compilers) written
8789
SIGNATURE OF TY COON, 1 April 1989
8790
Ty Coon, President of Vice
8792
This General Public License does not permit incorporating your program
8793
into proprietary programs. If your program is a subroutine library,
8794
you may consider it more useful to permit linking proprietary
8795
applications with the library. If this is what you want to do, use the
8796
GNU Library General Public License instead of this License.
8799
File: rheolef.info, Node: Concept Index, Up: Top