/[escript]/trunk/paso/src/Functions.c
ViewVC logotype

Annotation of /trunk/paso/src/Functions.c

Parent Directory Parent Directory | Revision Log Revision Log


Revision 1811 - (hide annotations)
Thu Sep 25 23:11:13 2008 UTC (11 years ago) by ksteube
File MIME type: text/plain
File size: 2154 byte(s)
Copyright updated in all files

1 gross 1476
2     /*******************************************************
3 ksteube 1811 *
4     * Copyright (c) 2003-2008 by University of Queensland
5     * Earth Systems Science Computational Center (ESSCC)
6     * http://www.uq.edu.au/esscc
7     *
8     * Primary Business: Queensland, Australia
9     * Licensed under the Open Software License version 3.0
10     * http://www.opensource.org/licenses/osl-3.0.php
11     *
12     *******************************************************/
13 gross 1476
14 ksteube 1811
15 gross 1476 #include "Common.h"
16     #include "Functions.h"
17 gross 1639 #include "PasoUtil.h"
18 phornby 1643 #include "Solver.h"
19 gross 1476 /*
20     * numerical calculation of the directional derivative J0w if F at x0 in the direction w. f0 is the value of F at x0.
21     * setoff is workspace
22     */
23    
24 gross 1804 err_t Paso_FunctionDerivative(double* J0w, const double* w, Paso_Function* F, const double *f0, const double *x0, double* setoff, const bool_t w_is_normalized, Paso_Performance* pp)
25 gross 1476 {
26     err_t err=0;
27 gross 1639 dim_t n=F->n;
28 gross 1476 double norm_w,epsnew,norm_x0;
29     epsnew=10.*sqrt(EPSILON);
30 gross 1798 if (! w_is_normalized) {
31     norm_w=Paso_l2(n,w,F->mpi_info);
32     } else {
33     norm_w=1.;
34     }
35 gross 1476
36     if (norm_w>0) {
37     epsnew = epsnew/norm_w;
38     norm_x0=Paso_l2(n,x0,F->mpi_info);
39     if (norm_x0>0) epsnew*=norm_x0;
40     Paso_LinearCombination(n,setoff,1.,x0,epsnew,w);
41 gross 1804 err=Paso_FunctionCall(F,J0w,setoff,pp);
42 gross 1639 if (err==NO_ERROR) {
43 gross 1476 Paso_Update(n,1./epsnew,J0w,-1./epsnew,f0); /* J0w = (J0w - f0)/epsnew; */
44     }
45 gross 1798 } else {
46     Paso_zeroes(n,J0w);
47 gross 1476 }
48     return err;
49     }
50    
51 gross 1639 /*
52     * sets value=F(arg)
53     *
54     */
55 gross 1804 err_t Paso_FunctionCall(Paso_Function * F,double* value, const double* arg, Paso_Performance *pp)
56 gross 1476 {
57 gross 1639 if (F!=NULL) {
58     switch(F->kind) {
59     case LINEAR_SYSTEM:
60 gross 1804 return Paso_Function_LinearSystem_call(F, value, arg,pp);
61 gross 1639 default:
62     return SYSTEM_ERROR;
63     }
64     }
65 phornby 1643 /* Added by PGH, assume a null pointe is an error */
66     return SYSTEM_ERROR;
67 gross 1639 }
68     /*
69     * clear Paso_Function
70     */
71     void Paso_Function_free(Paso_Function * F) {
72     if (F!=NULL) {
73 gross 1476
74 gross 1639 switch(F->kind) {
75     case LINEAR_SYSTEM:
76     Paso_Function_LinearSystem_free(F);
77     break;
78     default:
79     MEMFREE(F);
80     }
81     }
82 gross 1476 }

  ViewVC Help
Powered by ViewVC 1.1.26