On Thu, Dec 18, 2008 at 04:16:47PM +0100, Thomas Gazzola wrote: > Hello everyone! > > Since I post for the first time on the M-list, I thank you all for the > fenics suite. > > I am developing some code in C++, and I use dolfin in order to solve > some BVP problem. I did include dolfin.h and two forms compiled with > ffc in each file. I had a trouble when linking my stuff. > I am no C++ expert, so I don't know wether it is because of me or not. > > I did use the recent revision number 5349 (result of `hg id -n`). > > When compiling my stuff, I did have some _linking_ problems like: > > (...snip...) > Output.o: in function « > boost::numeric::ublas::compressed_matrix<double, > boost::numeric::ublas::basic_row_major<unsigned long, long>, 0ul, > boost::numeric::ublas::unbounded_array<unsigned long, > std::allocator<unsigned long> >, > boost::numeric::ublas::unbounded_array<double, std::allocator<double> > > >::size1() const»: > $PREFIX/local/include/dolfin/fem/assemble.h:89: multiple definitions > of « dolfin::assemble_system(dolfin::GenericMatrix&, > dolfin::GenericVector&, dolfin::Form const&, dolfin::Form const&, > std::vector<dolfin::DirichletBC const*, > std::allocator<dolfin::DirichletBC const*> >&, > dolfin::MeshFunction<unsigned int> const*, > dolfin::MeshFunction<unsigned int> const*, > dolfin::MeshFunction<unsigned int> const*, dolfin::GenericVector > const*, bool)» > main.o:/homes/tgax/local/include/dolfin/fem/assemble.h:89: first definition > here > (...snip...) > > In other words, I have troubles with function assemble_system in file > dolfin/fem/assemble.h . > > I had a look at the assemble.h file. It is a file with code calling > methods from the Assembler class (defined in files > dolfin/fem/Assembler.*). > > A `nm` on my *.o did show undefined symbols for function > assemble_system, in several object files. Assembler::assemble_system > were OK. > > The code from assemble.h seems to be compiled in several of my object > files, so the linker just cannot link. > > In order to solve this, I did split dolfin/fem/assemble.h into a .cpp > and a real header file. The SConstruct process automatically compiles > assemble.cpp, and the link could be done to my program. > > I don't know if it is a real problem of dolfin, or if I did something > wrong, but I wanted to let you known. What do you think? > > With the best regards, > > Thomas G.
> 13,15d12
> < #ifndef __ASSEMBLE_H
> < #define __ASSEMBLE_H
> <
> 27c24,27
> < bool reset_tensor);
> ---
> > bool reset_tensor=true)
> > {
> > Assembler::assemble(A, a, reset_tensor);
> > }
> 33c33,36
> < bool reset_tensor);
> ---
> > bool reset_tensor=true)
> > {
> > Assembler::assemble(A, a, sub_domain, reset_tensor);
> > }
> 41c44,47
> < bool reset_tensor);
> ---
> > bool reset_tensor=true)
> > {
> > Assembler::assemble(A, a, cell_domains, exterior_facet_domains,
> > interior_facet_domains);
> > }
> 49c55,58
> < bool reset_tensors);
> ---
> > bool reset_tensors=true)
> > {
> > Assembler::assemble_system(A, b, a, L, bc, reset_tensors);
> > }
> 57c66,69
> < bool reset_tensors);
> ---
> > bool reset_tensors=true)
> > {
> > Assembler::assemble_system(A, b, a, L, bcs, reset_tensors);
> > }
> 69c81,86
> < bool reset_tensors);
> ---
> > bool reset_tensors=true)
> > {
> > Assembler::assemble_system(A, b, a, L, bcs,
> > cell_domains, exterior_facet_domains,
> > interior_facet_domains,
> > x0, reset_tensors);
> > }
> 75c92,98
> < bool reset_tensor);
> ---
> > bool reset_tensor=true)
> > {
> > if (a.rank() != 0) error("Unable to assemble, form is not scalar.");
> > Scalar s;
> > Assembler::assemble(s, a, reset_tensor);
> > return s;
> > }
> 80c103,109
> < bool reset_tensor);
> ---
> > bool reset_tensor=true)
> > {
> > if (a.rank() != 0) error("Unable to assemble, form is not scalar.");
> > Scalar s;
> > Assembler::assemble(s, a, sub_domain, reset_tensor);
> > return s;
> > }
> 87,88c116,122
> < bool reset_tensor);
> < }
> ---
> > bool reset_tensor=true)
> > {
> > if (a.rank() != 0) error("Unable to assemble, form is not scalar.");
> > Scalar s;
> > Assembler::assemble(s, a, cell_domains, exterior_facet_domains,
> > interior_facet_domains, reset_tensor);
> > return s;
> > }
> 90c124
> < #endif
> ---
> > }
> // Copyright (C) 2007-2008 Anders Logg.
> // Licensed under the GNU LGPL Version 2.1.
> //
> // Modified by Garth N. Wells, 2008.
> //
> // First added: 2007-01-17
> // Last changed: 2008-11-16
> //
> // This file duplicates the Assembler::assemble* functions in
> // namespace dolfin, and adds special versions returning the value
> // directly for scalars. For documentation, refer to Assemble.h.
>
> #include <dolfin/log/log.h>
> #include <dolfin/la/Scalar.h>
> #include "Form.h"
> #include "Assembler.h"
> #include "assemble.h"
>
> namespace dolfin
> {
>
> using namespace dolfin;
> //--- Copies of assembly functions in Assembler.h ---
>
> /// Assemble tensor
> void assemble(GenericTensor& A,
> const Form& a,
> bool reset_tensor=true)
> {
> Assembler::assemble(A, a, reset_tensor);
> }
>
> /// Assemble tensor on sub domain
> void assemble(GenericTensor& A,
> const Form& a,
> const SubDomain& sub_domain,
> bool reset_tensor=true)
> {
> Assembler::assemble(A, a, sub_domain, reset_tensor);
> }
>
> /// Assemble tensor on sub domains
> void assemble(GenericTensor& A,
> const Form& a,
> const MeshFunction<uint>* cell_domains,
> const MeshFunction<uint>* exterior_facet_domains,
> const MeshFunction<uint>* interior_facet_domains,
> bool reset_tensor=true)
> {
> Assembler::assemble(A, a, cell_domains, exterior_facet_domains,
> interior_facet_domains);
> }
>
> /// Assemble system (A, b) and apply Dirichlet boundary condition
> void assemble_system(GenericMatrix& A,
> GenericVector& b,
> const Form& a,
> const Form& L,
> const DirichletBC& bc,
> bool reset_tensors=true)
> {
> Assembler::assemble_system(A, b, a, L, bc, reset_tensors);
> }
>
> /// Assemble system (A, b) and apply Dirichlet boundary conditions
> void assemble_system(GenericMatrix& A,
> GenericVector& b,
> const Form& a,
> const Form& L,
> std::vector<const DirichletBC*>& bcs,
> bool reset_tensors=true)
> {
> Assembler::assemble_system(A, b, a, L, bcs, reset_tensors);
> }
>
> /// Assemble system (A, b) on sub domains and apply Dirichlet boundary
> conditions
> void assemble_system(GenericMatrix& A,
> GenericVector& b,
> const Form& a,
> const Form& L,
> std::vector<const DirichletBC*>& bcs,
> const MeshFunction<uint>* cell_domains,
> const MeshFunction<uint>* exterior_facet_domains,
> const MeshFunction<uint>* interior_facet_domains,
> const GenericVector* x0,
> bool reset_tensors=true)
> {
> Assembler::assemble_system(A, b, a, L, bcs,
> cell_domains, exterior_facet_domains,
> interior_facet_domains,
> x0, reset_tensors);
> }
>
> //--- Specialized versions for scalars ---
>
> /// Assemble scalar
> double assemble(const Form& a,
> bool reset_tensor=true)
> {
> if (a.rank() != 0) error("Unable to assemble, form is not scalar.");
> Scalar s;
> Assembler::assemble(s, a, reset_tensor);
> return s;
> }
>
> /// Assemble scalar on sub domain
> double assemble(const Form& a,
> const SubDomain& sub_domain,
> bool reset_tensor=true)
> {
> if (a.rank() != 0) error("Unable to assemble, form is not scalar.");
> Scalar s;
> Assembler::assemble(s, a, sub_domain, reset_tensor);
> return s;
> }
>
> /// Assemble scalar on sub domains
> double assemble(const Form& a,
> const MeshFunction<uint>* cell_domains,
> const MeshFunction<uint>* exterior_facet_domains,
> const MeshFunction<uint>* interior_facet_domains,
> bool reset_tensor=true)
> {
> if (a.rank() != 0) error("Unable to assemble, form is not scalar.");
> Scalar s;
> Assembler::assemble(s, a, cell_domains, exterior_facet_domains,
> interior_facet_domains, reset_tensor);
> return s;
> }
>
> }
Sounds like a missing #ifndef in assemble.h.
I have fixed this now, see if it works.
--
Anders
signature.asc
Description: Digital signature
_______________________________________________ DOLFIN-dev mailing list [email protected] http://www.fenics.org/mailman/listinfo/dolfin-dev
