2
/// This file is part of Rheolef.
4
/// Copyright (C) 2000-2009 Pierre Saramito <Pierre.Saramito@imag.fr>
6
/// Rheolef is free software; you can redistribute it and/or modify
7
/// it under the terms of the GNU General Public License as published by
8
/// the Free Software Foundation; either version 2 of the License, or
9
/// (at your option) any later version.
11
/// Rheolef is distributed in the hope that it will be useful,
12
/// but WITHOUT ANY WARRANTY; without even the implied warranty of
13
/// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
14
/// GNU General Public License for more details.
16
/// You should have received a copy of the GNU General Public License
17
/// along with Rheolef; if not, write to the Free Software
18
/// Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
20
/// =========================================================================
22
// intel c++ portability issue:
24
// point.h declares friend norm(point x) member and this implementations
25
// generates a buggy result. It has been replaced by an extern implementation
26
// of norm(point x) and the result is portable.
28
#include "rheolef/point.h"
30
using namespace rheolef;
31
Float phi_1 (const point& x) { return norm(x) - 1; }
32
Float phi_2 (const point& x) { return sqrt(sqr(x[0])+sqr(x[1])+sqr(x[2])) - 1; }
33
int main (int argc, char**argv) {
34
#ifdef _RHEOLEF_USE_NEW_CODE
35
environment rheolef(argc, argv);
36
#endif // _RHEOLEF_USE_NEW_CODE
37
Float tol = sqrt(std::numeric_limits<Float>::epsilon());
39
Float p = sqrt(Float(14))-1;
42
Float err = fabs(p-p1);
44
cout << "phi (x)="<<p <<endl;
45
cout << "phi_2(x)="<<p2<<endl;
46
cout << "phi_1(x)="<<p1<<endl;
47
cout << "err ="<<err<<endl;
49
return (err < tol) ? 0 : 1;