tettri.c

Go to the documentation of this file.
00001 /*
00002     This file is part of magpar.
00003 
00004     Copyright (C) 2002-2010 Werner Scholz
00005 
00006     www:   http://www.magpar.net/
00007     email: magpar(at)magpar.net
00008 
00009     magpar is free software; you can redistribute it and/or modify
00010     it under the terms of the GNU General Public License as published by
00011     the Free Software Foundation; either version 2 of the License, or
00012     (at your option) any later version.
00013 
00014     magpar is distributed in the hope that it will be useful,
00015     but WITHOUT ANY WARRANTY; without even the implied warranty of
00016     MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
00017     GNU General Public License for more details.
00018 
00019     You should have received a copy of the GNU General Public License
00020     along with magpar; if not, write to the Free Software
00021     Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307  USA
00022 */
00023 
00024 static char const Id[] = "$Id: tettri.c 2962 2010-02-04 19:50:44Z scholz $\n\n";
00025 static char const Td[] = "$Today: " __FILE__ "  " __DATE__ "  "  __TIME__  " $\n\n";
00026 
00027 #include "util.h"
00028 
00029 PetscReal triarea(PetscReal *v0,PetscReal *v1,PetscReal *v2,PetscReal *n)
00030 {
00031   PetscReal a[ND],b[ND],outer[ND];
00032   PetscReal result;
00033 
00034   /* calculate a=v1-v0 */
00035   my_dcopy(ND,v1,1,a,1);
00036   my_daxpy(ND,-1.0,v0,1,a,1);
00037 
00038   /* calculate b=v2-v0 */
00039   my_dcopy(ND,v2,1,b,1);
00040   my_daxpy(ND,-1.0,v0,1,b,1);
00041 
00042   /* result=0.5*|a x b| */
00043   douter(ND,a,b,outer);
00044   result=0.5*my_dnrm2(ND,outer,1);
00045 
00046   if (n!=PETSC_NULL) {
00047     my_dscal(ND,1.0/(2.0*result),outer,1);
00048     my_dcopy(ND,outer,1,n,1);
00049   }
00050 
00051   return(result);
00052 }
00053 
00054 
00055 PetscReal tetvol(PetscReal *x1,PetscReal *x2,PetscReal *x3,PetscReal *x4)
00056 {
00057   PetscReal result;
00058 
00059   /* Calculate volume of the element */
00060   result=(1.0/6.0)*
00061     (x3[0]*x2[1]*x1[2] - x4[0]*x2[1]*x1[2] - x2[0]*x3[1]*x1[2] + x4[0]*x3[1]*x1[2] +
00062      x2[0]*x4[1]*x1[2] - x3[0]*x4[1]*x1[2] - x3[0]*x1[1]*x2[2] + x4[0]*x1[1]*x2[2] +
00063      x1[0]*x3[1]*x2[2] - x4[0]*x3[1]*x2[2] - x1[0]*x4[1]*x2[2] + x3[0]*x4[1]*x2[2] +
00064      x2[0]*x1[1]*x3[2] - x4[0]*x1[1]*x3[2] - x1[0]*x2[1]*x3[2] + x4[0]*x2[1]*x3[2] +
00065      x1[0]*x4[1]*x3[2] - x2[0]*x4[1]*x3[2] - x2[0]*x1[1]*x4[2] + x3[0]*x1[1]*x4[2] +
00066      x1[0]*x2[1]*x4[2] - x3[0]*x2[1]*x4[2] - x1[0]*x3[1]*x4[2] + x2[0]*x3[1]*x4[2]
00067     );
00068 
00069   return(result);
00070 }
00071 
00072 int tetnb(int *tet1, int *tet2, int *vshared)
00073 {
00074   int t1[NV],t2[NV];
00075 
00076   for (int i=0;i<NV;i++) {
00077     t1[i]=tet1[i];
00078     t2[i]=tet2[i];
00079   }
00080 
00081   /* sort nodes */
00082   PetscSortInt(NV,(PetscInt*)t1);
00083   PetscSortInt(NV,(PetscInt*)t2);
00084 
00085   int cntshared=0;
00086 
00087   for (int i=0;i<NV;i++) {
00088     for (int j=0;j<NV;j++) {
00089       if (t1[i]==t2[j]) {
00090         if (vshared!=PETSC_NULL) vshared[cntshared]=t1[i];
00091         cntshared++;
00092       }
00093     }
00094   }
00095 
00096   assert(cntshared<=NV);
00097 
00098   return(cntshared);
00099 }
00100 
00101 int TetSharedVertices(int *v1,int *v2,int *s)
00102 {
00103   int k;
00104   k=0;
00105   for (int i=0;i<NV;i++) {
00106     for (int j=0;j<NV;j++) {
00107       if (v1[i]==v2[j]) {
00108         s[k]=v1[i];
00109         k++;
00110         continue;
00111       }
00112     }
00113   }
00114   return(k);
00115 }

magpar - Parallel Finite Element Micromagnetics Package
Copyright (C) 2002-2010 Werner Scholz