# Contents of /trunk/paso/src/Functions.c

Revision 1804 - (show annotations)
Wed Sep 24 07:52:19 2008 UTC (11 years ago) by gross
File MIME type: text/plain
File size: 2141 byte(s)
```a robister version of the upwinding scheme
```
 1 /* \$Id:\$ */ 2 3 /******************************************************* 4 * 5 * Copyright 2008 by University of Queensland 6 * 7 * http://esscc.uq.edu.au 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 14 #include "Common.h" 15 #include "Functions.h" 16 #include "PasoUtil.h" 17 #include "Solver.h" 18 /* 19 * numerical calculation of the directional derivative J0w if F at x0 in the direction w. f0 is the value of F at x0. 20 * setoff is workspace 21 */ 22 23 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) 24 { 25 err_t err=0; 26 dim_t n=F->n; 27 double norm_w,epsnew,norm_x0; 28 epsnew=10.*sqrt(EPSILON); 29 if (! w_is_normalized) { 30 norm_w=Paso_l2(n,w,F->mpi_info); 31 } else { 32 norm_w=1.; 33 } 34 35 if (norm_w>0) { 36 epsnew = epsnew/norm_w; 37 norm_x0=Paso_l2(n,x0,F->mpi_info); 38 if (norm_x0>0) epsnew*=norm_x0; 39 Paso_LinearCombination(n,setoff,1.,x0,epsnew,w); 40 err=Paso_FunctionCall(F,J0w,setoff,pp); 41 if (err==NO_ERROR) { 42 Paso_Update(n,1./epsnew,J0w,-1./epsnew,f0); /* J0w = (J0w - f0)/epsnew; */ 43 } 44 } else { 45 Paso_zeroes(n,J0w); 46 } 47 return err; 48 } 49 50 /* 51 * sets value=F(arg) 52 * 53 */ 54 err_t Paso_FunctionCall(Paso_Function * F,double* value, const double* arg, Paso_Performance *pp) 55 { 56 if (F!=NULL) { 57 switch(F->kind) { 58 case LINEAR_SYSTEM: 59 return Paso_Function_LinearSystem_call(F, value, arg,pp); 60 default: 61 return SYSTEM_ERROR; 62 } 63 } 64 /* Added by PGH, assume a null pointe is an error */ 65 return SYSTEM_ERROR; 66 } 67 /* 68 * clear Paso_Function 69 */ 70 void Paso_Function_free(Paso_Function * F) { 71 if (F!=NULL) { 72 73 switch(F->kind) { 74 case LINEAR_SYSTEM: 75 Paso_Function_LinearSystem_free(F); 76 break; 77 default: 78 MEMFREE(F); 79 } 80 } 81 }