Mercurial > fem-fenics-eugenio

comparison src/feval.cc @ 241:0f14cdbcaed3

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
feval.cc: stop computation after an error occurs
author Eugenio Gianniti <eugenio.gianniti@mail.polimi.it>
date Wed, 25 Jun 2014 20:36:04 +0200
parents 418a5119047b
children 5e9b5bbdc56b
comparison
equal deleted inserted replaced
240:a57001686abb 241:0f14cdbcaed3
1 /* 1 /*
2 Copyright (C) 2013 Marco Vassallo <gedeone-octave@users.sourceforge.net> 2 Copyright (C) 2013 Marco Vassallo <gedeone-octave@users.sourceforge.net>
3 Copyright (C) 2014 Eugenio Gianniti <eugenio.gianniti@mail.polimi.it>
3 4
4 This program is free software; you can redistribute it and/or modify it under 5 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 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 Foundation; either version 3 of the License, or (at your option) any later
7 version. 8 version.
78 if (dims != newdims) 79 if (dims != newdims)
79 { 80 {
80 std::string msg = "all the input matrices should"; 81 std::string msg = "all the input matrices should";
81 msg += " have the same size"; 82 msg += " have the same size";
82 error (msg.c_str ()); 83 error (msg.c_str ());
84 break;
83 } 85 }
84 } 86 }
85 coordinates.push_back (aux); 87 coordinates.push_back (aux);
86 } 88 }
87 } 89 }
88 90
89 octave_idx_type vdim = f->value_dimension (0); 91 if (! error_state)
90 if (nargout != vdim) 92 {
91 error ("wrong number of output arguments"); 93 octave_idx_type vdim = f->value_dimension (0);
94 if (nargout != vdim)
95 error ("wrong number of output arguments");
96 else
97 {
98 std::vector <Matrix> evaluations;
99 for (octave_idx_type out = 0; out < vdim; ++out)
100 evaluations.push_back (Matrix (dims));
92 101
93 std::vector <Matrix> evaluations; 102 for (octave_idx_type k = 0; k < dims.numel (); ++k)
94 for (octave_idx_type out = 0; out < vdim; ++out) 103 {
95 evaluations.push_back (Matrix (dims)); 104 Array<double> point (dim_vector (pdim, 1));
105 for (octave_idx_type el = 0; el < pdim; ++el)
106 point (el) = coordinates[el] (k);
107 dolfin::Array<double>
108 x (point.length (), point.fortran_vec ());
96 109
97 for (octave_idx_type k = 0; k < dims.numel (); ++k) 110 Array<double> res (dim_vector (vdim, 1));
98 { 111 dolfin::Array<double>
99 Array<double> point (dim_vector (pdim, 1)); 112 values (res.length (), res.fortran_vec ());
100 for (octave_idx_type el = 0; el < pdim; ++el) 113 try
101 point (el) = coordinates[el] (k); 114 {
102 dolfin::Array<double> 115 f->eval (values, x);
103 x (point.length (), point.fortran_vec ()); 116 }
117 catch (std::runtime_error & err)
118 {
119 std::string msg = "cannot evaluate a function";
120 msg += " outside of its domain";
121 error (msg.c_str ());
122 break;
123 }
104 124
105 Array<double> res (dim_vector (vdim, 1)); 125 for (octave_idx_type el = 0; el < vdim; ++el)
106 dolfin::Array<double> 126 evaluations[el] (k) = res (el);
107 values (res.length (), res.fortran_vec ()); 127 }
108 try 128
109 { 129 if (! error_state)
110 f->eval (values, x); 130 {
131 for (std::vector<Matrix>::iterator it =
132 evaluations.begin ();
133 it != evaluations.end (); ++it)
134 retval.append (octave_value (*it));
135 }
111 } 136 }
112 catch (std::runtime_error & err)
113 {
114 std::string msg = "cannot evaluate a function";
115 msg += " outside of its domain";
116 error (msg.c_str ());
117 }
118
119 for (octave_idx_type el = 0; el < vdim; ++el)
120 evaluations[el] (k) = res (el);
121 } 137 }
122
123 for (std::vector<Matrix>::iterator it = evaluations.begin ();
124 it != evaluations.end (); ++it)
125 retval.append (octave_value (*it));
126 } 138 }
127 } 139 }
128 } 140 }
129 141
130 return retval; 142 return retval;