00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017
00018
00019
00020
00021
00022
00023
00024
00025
00026
00027
00028 #include "util.h"
00029 #include "petscblaslapack.h"
00030
00031
00032
00033
00034
00035 #if !defined(PETSC_USE_COMPLEX)
00036
00037 #if defined(PETSC_USES_FORTRAN_SINGLE)
00038 #define DGESV SGESV
00039 #endif
00040
00041 #if defined(PETSC_USE_SINGLE)
00042
00043 #if defined(PETSC_HAVE_FORTRAN_UNDERSCORE) || defined(PETSC_BLASLAPACK_F2C)
00044 #define LAgesv_ sgesv_
00045 #elif defined(PETSC_HAVE_FORTRAN_CAPS)
00046 #define LAgesv_ SGESV
00047 #else
00048 #define LAgesv_ sgesv
00049 #endif
00050
00051 #else
00052
00053 #if defined(PETSC_HAVE_FORTRAN_UNDERSCORE) || defined(PETSC_BLASLAPACK_F2C)
00054 #define LAgesv_ dgesv_
00055 #elif defined(PETSC_HAVE_FORTRAN_CAPS)
00056 #define LAgesv_ DGESV
00057 #else
00058 #define LAgesv_ dgesv
00059 #endif
00060
00061 #endif
00062
00063 #else
00064
00065 #if defined(PETSC_USES_FORTRAN_SINGLE)
00066 #define ZGESV CGESV
00067 #endif
00068
00069 #if defined(PETSC_HAVE_FORTRAN_UNDERSCORE) || defined(PETSC_BLASLAPACK_F2C)
00070 #define LAgesv_ zgesv_
00071 #elif defined(PETSC_HAVE_FORTRAN_CAPS)
00072 #define LAgesv_ ZGESV
00073 #else
00074 #define LAgesv_ zgesv
00075 #endif
00076
00077 #endif
00078
00079 EXTERN_C_BEGIN
00080 EXTERN void LAgesv_(const PetscBLASInt*, const PetscBLASInt*, PetscScalar*, const PetscBLASInt*, PetscBLASInt*, PetscScalar*, const PetscBLASInt*, PetscBLASInt*);
00081 EXTERN_C_END
00082
00083
00084
00085
00086
00087
00088
00089
00090 int barycent(PetscReal *p,PetscReal *a,PetscReal *b,PetscReal *c,PetscReal *d,PetscReal *bary)
00091 {
00092 my_dcopy(ND,p,1,bary,1);
00093 my_daxpy(ND,-1.0,d,1,bary,1);
00094
00095
00096 PetscReal tetmat[ND*ND];
00097 my_dcopy(ND,a,1,tetmat,1); my_daxpy(ND,-1.0,d,1,tetmat,1);
00098 my_dcopy(ND,b,1,tetmat+ND,1); my_daxpy(ND,-1.0,d,1,tetmat+ND,1);
00099 my_dcopy(ND,c,1,tetmat+2*ND,1);my_daxpy(ND,-1.0,d,1,tetmat+2*ND,1);
00100
00101
00102 int one=1;
00103 int three=ND;
00104 int ipiv[ND];
00105 int info;
00106 LAgesv_(&three,&one,tetmat,&three,ipiv,bary,&three,&info);CHKERRQ(info);
00107
00108 if (bary[0]>=0.0 && bary[1]>=0.0 && bary[2]>=0.0 &&
00109 (bary[0]+bary[1]+bary[2])<1.0) {
00110 return(1);
00111 }
00112 else {
00113 return(0);
00114 }
00115 }
00116