Mercurial > fem-fenics-eugenio
comparison src/mark.cc @ 267:53039ac90368
Mark meshfunction using subdomain
* src/SubDomain.cc: return subdomain to use for marking
* src/MeshFunction.cc: return a user defined meshfunction
* src/mark.cc: with the information from subdomain, mark a meshfunction
author | Eugenio Gianniti <eugenio.gianniti@mail.polimi.it> |
---|---|
date | Tue, 12 Aug 2014 15:42:50 +0200 |
parents | |
children | 61830a4f9ab9 |
comparison
equal
deleted
inserted
replaced
266:6b37560b7cbb | 267:53039ac90368 |
---|---|
1 /* | |
2 Copyright (C) 2014 Eugenio Gianniti <eugenio.gianniti@mail.polimi.it> | |
3 | |
4 This program is free software; you can redistribute it and/or modify it under | |
5 the terms of the GNU General Public License as published by the Free Software | |
6 Foundation; either version 3 of the License, or (at your option) any later | |
7 version. | |
8 | |
9 This program is distributed in the hope that it will be useful, but WITHOUT | |
10 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or | |
11 FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License for more | |
12 details. | |
13 | |
14 You should have received a copy of the GNU General Public License along with | |
15 this program; if not, see <http://www.gnu.org/licenses/>. | |
16 */ | |
17 | |
18 #include <dolfin.h> | |
19 #include "meshfunction.h" | |
20 #include "subdomain.h" | |
21 #include "dolfin_compat.h" | |
22 | |
23 DEFUN_DLD (mark, args, nargout, | |
24 "-*- texinfo -*-\n\ | |
25 @deftypefn {Function File} {[@var{markers}]} = \ | |
26 mark (@var{subdomain}, @var{meshfunction}, @var{number})\n\ | |
27 Mark @var{meshfunction} with @var{number} on the entities contained in \ | |
28 @var{subdomain}. The output, @var{markers}, is a copy of @var{meshfunction} \ | |
29 where entities inside @var{subdomain} are marked.\n\ | |
30 @seealso{MeshFunction, SubDomain}\n\ | |
31 @end deftypefn") | |
32 { | |
33 int nargin = args.length (); | |
34 octave_value retval; | |
35 | |
36 if (nargin < 3 || nargin > 3 ||nargout > 1) | |
37 print_usage (); | |
38 else | |
39 { | |
40 if (! meshfunction_type_loaded) | |
41 { | |
42 meshfunction::register_type (); | |
43 meshfunction_type_loaded = true; | |
44 mlock (); | |
45 } | |
46 if (! subdomain_type_loaded) | |
47 { | |
48 subdomain::register_type (); | |
49 subdomain_type_loaded = true; | |
50 mlock (); | |
51 } | |
52 | |
53 if (args(0).type_id () == subdomain::static_type_id () && | |
54 args(1).type_id () == meshfunction::static_type_id () && | |
55 args(2).is_real_scalar ()) | |
56 { | |
57 std::size_t number = args(2).ulong_value (); | |
58 subdomain const & sd = static_cast <subdomain const &> | |
59 (args(0).get_rep ()); | |
60 meshfunction const & mf = static_cast <meshfunction const &> | |
61 (args(1).get_rep ()); | |
62 | |
63 if (! error_state) | |
64 { | |
65 dolfin::MeshFunction <std::size_t> out = mf.get_mf (); | |
66 SHARED_PTR <subdomain_rep const> const & sdrep = sd.get_psd (); | |
67 | |
68 sdrep->mark (out, number); | |
69 retval = new meshfunction (mf.get_str (), out); | |
70 } | |
71 } | |
72 else | |
73 error ("invalid input arguments"); | |
74 } | |
75 | |
76 return retval; | |
77 } |