1c6db04a5SJed Brown#include <petscsys.h> 2c4da4110SBarry Smith 3*887451b3SStefano Zampini#define NOTLIKELY PETSC_INFINITY 4c4da4110SBarry Smith 5*887451b3SStefano ZampiniPETSC_EXTERN PetscErrorCode PFApply_String(void *value,PetscInt n,const PetscScalar *in,PetscScalar *out) 6c4da4110SBarry Smith{ 7f236b2adSBarry Smith PetscInt i; 8f236b2adSBarry Smith PetscScalar x,y,z,f = NOTLIKELY,x1 = 0,x2 = 0,x3 = 0,x4 = 0,x5 = 0; 9f236b2adSBarry Smith PetscScalar f1 = NOTLIKELY, f2 = 0,f3 = 0,f4 = 0,f5 = 0; 10c4da4110SBarry Smith 11f236b2adSBarry Smith (void)x; 12f236b2adSBarry Smith (void)y; 13f236b2adSBarry Smith (void)z; 14f236b2adSBarry Smith (void)f; 15f236b2adSBarry Smith (void)x1; 16f236b2adSBarry Smith (void)x2; 17f236b2adSBarry Smith (void)x3; 18f236b2adSBarry Smith (void)x4; 19f236b2adSBarry Smith (void)x5; 20f236b2adSBarry Smith (void)f1; 21f236b2adSBarry Smith (void)f2; 22f236b2adSBarry Smith (void)f3; 23f236b2adSBarry Smith (void)f4; 24f236b2adSBarry Smith (void)f5; 25c4da4110SBarry Smith for (i=0; i<n; i++) { 26f236b2adSBarry Smith x1 = x = in[_NIN_*i]; 27f236b2adSBarry Smith#if (_NIN_ > 1) 28f236b2adSBarry Smith x2 = y = in[_NIN_*i+1]; 29c4da4110SBarry Smith#endif 30f236b2adSBarry Smith#if (_NIN_ > 2) 31f236b2adSBarry Smith x3 = z = in[_NIN_*i+2]; 32c4da4110SBarry Smith#endif 33f236b2adSBarry Smith#if (_NIN_ > 3) 34f236b2adSBarry Smith x4 = in[_NIN_*i+3]; 35c4da4110SBarry Smith#endif 36f236b2adSBarry Smith#if (_NIN_ > 4) 37f236b2adSBarry Smith x5 = in[_NIN_*i+4]; 38c4da4110SBarry Smith#endif 39c4da4110SBarry Smith FUNCTION; 40c4da4110SBarry Smith if (f == NOTLIKELY) { 41c4da4110SBarry Smith out[_NOUT_*i] = f1; 42c4da4110SBarry Smith } else { 43c4da4110SBarry Smith out[_NOUT_*i] = f; 44c4da4110SBarry Smith } 45c4da4110SBarry Smith#if (_NOUT_ > 1) 46c4da4110SBarry Smith out[_NOUT_*i+1] = f2; 47c4da4110SBarry Smith#endif 48c4da4110SBarry Smith#if (_NOUT_ > 2) 49c4da4110SBarry Smith out[_NOUT_*i+2] = f3; 50c4da4110SBarry Smith#endif 51c4da4110SBarry Smith#if (_NOUT_ > 3) 52c4da4110SBarry Smith out[_NOUT_*i+3] = f4; 53c4da4110SBarry Smith#endif 54c4da4110SBarry Smith#if (_NOUT_ > 4) 55c4da4110SBarry Smith out[_NOUT_*i+4] = f5; 56c4da4110SBarry Smith#endif 57c4da4110SBarry Smith } 58*887451b3SStefano Zampini return(PETSC_SUCCESS); 59c4da4110SBarry Smith} 60