Mercurial > fem-fenics-eugenio

comparison src/feval.cc @ 268:61830a4f9ab9

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
Improve formatting
author Eugenio Gianniti <eugenio.gianniti@mail.polimi.it>
date Thu, 14 Aug 2014 12:26:55 +0200
parents 5e9b5bbdc56b
children
comparison
equal deleted inserted replaced
267:53039ac90368 268:61830a4f9ab9
37 @end deftypefn") 37 @end deftypefn")
38 { 38 {
39 39
40 int nargin = args.length (); 40 int nargin = args.length ();
41 octave_value_list retval; 41 octave_value_list retval;
42 42
43 if (! function_type_loaded) 43 if (! function_type_loaded)
44 { 44 {
45 function::register_type (); 45 function::register_type ();
46 function_type_loaded = true; 46 function_type_loaded = true;
47 mlock (); 47 mlock ();
48 } 48 }
49 49
50 if (args(0).type_id () == function::static_type_id ()) 50 if (args(0).type_id () == function::static_type_id ())
51 { 51 {
52 const function & fspo = 52 const function & fspo =
53 static_cast<const function&> (args(0).get_rep ()); 53 static_cast<const function &> (args(0).get_rep ());
54 54
55 if (!error_state) 55 if (!error_state)
56 { 56 {
57 const SHARED_PTR <const dolfin::Function> 57 const SHARED_PTR <const dolfin::Function>
58 & f = fspo.get_pfun (); 58 & f = fspo.get_pfun ();
59 59
60 octave_idx_type pdim = f->geometric_dimension (); 60 octave_idx_type pdim = f->geometric_dimension ();
61 if (nargin != pdim + 1) 61 if (nargin != pdim + 1)
62 print_usage (); 62 { print_usage (); }
63 else 63 else
64 { 64 {
65 std::vector <Matrix> coordinates; 65 std::vector <Matrix> coordinates;
66 dim_vector dims; 66 dim_vector dims;
67 67
68 for (octave_idx_type in = 1; in <= pdim; ++in) 68 for (octave_idx_type in = 1; in <= pdim; ++in)
69 { 69 {
70 if (! args(in).is_real_type ()) 70 if (! args(in).is_real_type ())
71 error ("invalid argument"); 71 { error ("invalid argument"); }
72 else 72 else
73 { 73 {
74 Matrix aux = args(in).matrix_value (); 74 Matrix aux = args(in).matrix_value ();
75 if (in == 1) 75 if (in == 1)
76 dims = aux.dims (); 76 { dims = aux.dims (); }
77 else 77 else
78 { 78 {
79 dim_vector newdims = aux.dims (); 79 dim_vector newdims = aux.dims ();
80 if (dims != newdims) 80 if (dims != newdims)
81 { 81 {
91 91
92 if (! error_state) 92 if (! error_state)
93 { 93 {
94 octave_idx_type vdim = f->value_dimension (0); 94 octave_idx_type vdim = f->value_dimension (0);
95 if (nargout != vdim) 95 if (nargout != vdim)
96 error ("wrong number of output arguments"); 96 { error ("wrong number of output arguments"); }
97 else 97 else
98 { 98 {
99 std::vector <Matrix> evaluations; 99 std::vector <Matrix> evaluations;
100 for (octave_idx_type out = 0; out < vdim; ++out) 100 for (octave_idx_type out = 0; out < vdim; ++out)
101 evaluations.push_back (Matrix (dims)); 101 { evaluations.push_back (Matrix (dims)); }
102 102
103 for (octave_idx_type k = 0; k < dims.numel (); ++k) 103 for (octave_idx_type k = 0; k < dims.numel (); ++k)
104 { 104 {
105 Array<double> point (dim_vector (pdim, 1)); 105 Array<double> point (dim_vector (pdim, 1));
106 for (octave_idx_type el = 0; el < pdim; ++el) 106 for (octave_idx_type el = 0; el < pdim; ++el)
107 point (el) = coordinates[el] (k); 107 { point (el) = coordinates[el] (k); }
108 dolfin::Array<double> 108 dolfin::Array<double>
109 x (point.length (), point.fortran_vec ()); 109 x (point.length (), point.fortran_vec ());
110 110
111 Array<double> res (dim_vector (vdim, 1)); 111 Array<double> res (dim_vector (vdim, 1));
112 dolfin::Array<double> 112 dolfin::Array<double>
113 values (res.length (), res.fortran_vec ()); 113 values (res.length (), res.fortran_vec ());
114 try 114 try
115 { 115 {
116 f->eval (values, x); 116 f->eval (values, x);
117 } 117 }
122 error (msg.c_str ()); 122 error (msg.c_str ());
123 break; 123 break;
124 } 124 }
125 125
126 for (octave_idx_type el = 0; el < vdim; ++el) 126 for (octave_idx_type el = 0; el < vdim; ++el)
127 evaluations[el] (k) = res (el); 127 { evaluations[el] (k) = res (el); }
128 } 128 }
129 129
130 if (! error_state) 130 if (! error_state)
131 { 131 {
132 for (std::vector<Matrix>::iterator it = 132 for (std::vector<Matrix>::iterator it =
133 evaluations.begin (); 133 evaluations.begin ();
134 it != evaluations.end (); ++it) 134 it != evaluations.end (); ++it)
135 retval.append (octave_value (*it)); 135 { retval.append (octave_value (*it)); }
136 } 136 }
137 } 137 }
138 } 138 }
139 } 139 }
140 } 140 }