1 /* 2 * Copyright (C) 1998, 2000-2007, 2010, 2011, 2012, 2013 SINTEF ICT, 3 * Applied Mathematics, Norway. 4 * 5 * Contact information: E-mail: tor.dokken@sintef.no 6 * SINTEF ICT, Department of Applied Mathematics, 7 * P.O. Box 124 Blindern, 8 * 0314 Oslo, Norway. 9 * 10 * This file is part of SISL. 11 * 12 * SISL is free software: you can redistribute it and/or modify 13 * it under the terms of the GNU Affero General Public License as 14 * published by the Free Software Foundation, either version 3 of the 15 * License, or (at your option) any later version. 16 * 17 * SISL is distributed in the hope that it will be useful, 18 * but WITHOUT ANY WARRANTY; without even the implied warranty of 19 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the 20 * GNU Affero General Public License for more details. 21 * 22 * You should have received a copy of the GNU Affero General Public 23 * License along with SISL. If not, see 24 * <http://www.gnu.org/licenses/>. 25 * 26 * In accordance with Section 7(b) of the GNU Affero General Public 27 * License, a covered work must retain the producer line in every data 28 * file that is created or manipulated using SISL. 29 * 30 * Other Usage 31 * You can be released from the requirements of the license by purchasing 32 * a commercial license. Buying such a license is mandatory as soon as you 33 * develop commercial activities involving the SISL library without 34 * disclosing the source code of your own applications. 35 * 36 * This file may be used in accordance with the terms contained in a 37 * written agreement between you and SINTEF ICT. 38 */ 39 40 #include "sisl-copyright.h" 41 42 /* 43 * 44 * $Id: sh1464.c,v 1.1 1994-04-21 12:10:42 boh Exp $ 45 * 46 */ 47 48 49 #define SH1464 50 51 #include "sislP.h" 52 53 #if defined(SISLNEEDPROTOTYPES) 54 typedef void (*fshapeProc)(double [],double [],int,int,int *); 55 #else 56 typedef void (*fshapeProc)(); 57 #endif 58 59 #if defined(SISLNEEDPROTOTYPES) 60 void sh1464(fshapeProc fshape,SISLCurve * vboundc[],int icurv,double etwist[],double etang[],double eder[],int * jstat)61 sh1464(fshapeProc fshape,SISLCurve *vboundc[],int icurv, 62 double etwist[],double etang[],double eder[],int *jstat) 63 #else 64 void sh1464(fshape,vboundc,icurv,etwist,etang,eder,jstat) 65 fshapeProc fshape; 66 double etwist[],etang[],eder[]; 67 SISLCurve *vboundc[]; 68 int icurv,*jstat; 69 #endif 70 /* 71 ********************************************************************* 72 * 73 * PURPOSE : Given a five sided vertex region, evaluate the first 74 * blending surface in the corner lying in the middle of the 75 * vertex region. Compute the tangent vectors in the middle 76 * vertex along the inner boundaries of the region. 77 * 78 * 79 * 80 * INPUT : fshape - Application driven routine that gives the user an 81 * ability to change the middle point of the region 82 * (the vertex at which the blending surfaces meet), 83 * and the tangent vectors in the middle point along 84 * the curves which divedes the region. 85 * vboundc - Position and cross-tangent curves around the vertex 86 * region. For each edge of the region position and cross- 87 * tangent curves are given. The curves follow each other 88 * around the region and are oriented counter-clock-wise. 89 * The dimension of the array is 10. 90 * icurv - Number of sides. icurv = 5. 91 * etwist - Twist-vectors of the corners of the vertex region. The 92 * first element of the array is the twist in the corner 93 * before the first edge, etc. The dimension of the array 94 * is icurve*kdim. 95 * 96 * 97 * OUTPUT : etang - Tangent vectors at the midpoint of the vertex region. 98 * The dimension is icurv*idim. 99 * eder - Value, first and second derivative of the first blending 100 * surface in the corner at the midpoint. The sequence is the 101 * following : Value, 1. derivative in 1. parameter direction, 102 * 1. derivative in the 2. parameter direction, 2. derivative 103 * in the 1. parameter direction, mixed derivative and 2. 104 * derivative in the 2. parameter direction. Dimension 6*idim. 105 * jstat - status messages 106 * > 0 : warning 107 * = 0 : ok 108 * < 0 : error 109 * 110 * 111 * METHOD : Evaluate the Gregory Charrot function in the midpoint of the 112 * vertex region. Compute the wanted derivatives. 113 * 114 * REFERENCES : 115 * 116 * USE : 3D geometry only. 117 * 118 *- 119 * CALLS : sh1467 - Evaluate Gregory Charrot function over 120 * 5-sided region. 121 * 122 * WRITTEN BY : Vibeke Skytt, SI, 03.90. 123 * 124 ********************************************************************* 125 */ 126 { 127 int kstat = 0; /* Status variable. */ 128 int kder = 2; /* Number of derviatives to evaluate. */ 129 int ki; /* Counter. */ 130 int kdim = 3; /* Dimension of geometry space. */ 131 double tlambda = (double)1.0/sqrt((double)5.0); /* Constant. */ 132 double tl1 = (double)2.0*tlambda*tan(PI/(double)5.0); /* Constant. */ 133 double tl2 = sin((double)0.3*PI); /* Constant. */ 134 double tconst1 = tl1/(double)2.0 - tlambda; /* Constant. */ 135 double tconst2 = tl2 - tlambda; /* Constant. */ 136 double sbar[5]; /* Barycentric coordinates of the blending function. */ 137 double sder[18]; /* Value and derivatives of blending function. */ 138 139 /* Set up the barycentric coordinates of the midpoint of the region. */ 140 141 sbar[0] = sbar[1] = sbar[2] = sbar[3] = sbar[4] = tlambda; 142 143 /* Evaluate the Gregory Charrot function at the midpoint. */ 144 145 sh1467(vboundc,etwist,kder,sbar,sder,&kstat); 146 if (kstat < 0) goto error; 147 148 /* Compute tangent vectors. */ 149 150 for (ki=0; ki<kdim; ki++) 151 { 152 etang[ki] = sder[kdim+ki]*tconst1 + sder[2*kdim+ki]*tconst2; 153 etang[kdim+ki] = -sder[kdim+ki]*tlambda + sder[2*kdim+ki]*tconst1; 154 etang[2*kdim+ki] = sder[kdim+ki]*tconst1 - sder[2*kdim+ki]*tlambda; 155 etang[3*kdim+ki] = sder[kdim+ki]*tconst2 + sder[2*kdim+ki]*tconst1; 156 etang[4*kdim+ki] = sder[kdim+ki]*tconst2 + sder[2*kdim+ki]*tconst2; 157 } 158 159 /* Application driven routine to alter the midpoint and tangents in the 160 midpoint. */ 161 162 fshape(sder,etang,kdim,icurv,&kstat); 163 if (kstat < 0) goto error; 164 165 /* Copy value and 1. derivatives of first patch. */ 166 167 memcopy(eder,sder,kdim,DOUBLE); 168 memcopy(eder+kdim,etang+4*kdim,kdim,DOUBLE); 169 memcopy(eder+2*kdim,etang,kdim,DOUBLE); 170 171 /* Compute 2. derivatives. */ 172 173 for (ki=0; ki<kdim; ki++) 174 { 175 eder[3*kdim+ki] = sder[3*kdim+ki]*tconst2*tconst2 176 + (double)2.0*sder[4*kdim+ki]*tconst2*tconst2 177 + sder[5*kdim+ki]*tconst2*tconst2; 178 eder[4*kdim+ki] = sder[3*kdim+ki]*tconst1*tconst2 179 + sder[4*kdim+ki]*tconst2*(tconst1+tconst2) 180 + sder[5*kdim+ki]*tconst2*tconst2; 181 eder[5*kdim+ki] = sder[3*kdim+ki]*tconst1*tconst1 182 - (double)2.0*sder[4*kdim+ki]*tconst1*tconst2 183 + sder[5*kdim+ki]*tconst2*tconst2; 184 } 185 186 *jstat = 0; 187 goto out; 188 189 /* Error in a lower level function. */ 190 191 error: 192 *jstat = kstat; 193 goto out; 194 195 out : 196 return; 197 } 198