1 /* -*- coding: utf-8 -*- */
2 /**
3 \ingroup anespec gshhs mate legendre
4 @{
5 \file mate.c
6 \brief Definición de funciones para la realización de cálculos matemáticos
7 generales.
8 \author José Luis García Pallero, jgpallero@gmail.com
9 \date 17 de mayo de 2010
10 \copyright
11 Copyright (c) 2009-2016, José Luis García Pallero. All rights reserved.
12 \par
13 Redistribution and use in source and binary forms, with or without modification,
14 are permitted provided that the following conditions are met:
15 \par
16 - Redistributions of source code must retain the above copyright notice, this
17 list of conditions and the following disclaimer.
18 - Redistributions in binary form must reproduce the above copyright notice, this
19 list of conditions and the following disclaimer in the documentation and/or
20 other materials provided with the distribution.
21 - Neither the name of the copyright holders nor the names of its contributors
22 may be used to endorse or promote products derived from this software without
23 specific prior written permission.
24 \par
25 THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS" AND
26 ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED
27 WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE
28 DISCLAIMED. IN NO EVENT SHALL COPYRIGHT HOLDER BE LIABLE FOR ANY DIRECT,
29 INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING,
30 BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE,
31 DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
32 LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE
33 OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF
34 ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
35 */
36 /******************************************************************************/
37 /******************************************************************************/
38 #include"libgeoc/mate.h"
39 /******************************************************************************/
40 /******************************************************************************/
41 /**
42 \brief Array que almacena de manera explícita los factoriales de los números que
43 van de \em 0 a #GEOC_MATE_CONST_LDBL_NMAXFAC.
44
45 Este array se basa en el array \p fact_table, que se puede encontrar en
46 el fichero \p gamma.c, de la biblioteca GSL.
47 */
48 static __mateFactExpl __tablaFacExpl[GEOC_MATE_CONST_LDBL_NMAXFAC+1]=
49 {
50 { 0, 1.0 },
51 { 1, 1.0 },
52 { 2, 2.0 },
53 { 3, 6.0 },
54 { 4, 24.0 },
55 { 5, 120.0 },
56 { 6, 720.0 },
57 { 7, 5040.0 },
58 { 8, 40320.0 },
59 { 9, 362880.0 },
60 { 10, 3628800.0 },
61 { 11, 39916800.0 },
62 { 12, 479001600.0 },
63 { 13, 6227020800.0 },
64 { 14, 87178291200.0 },
65 { 15, 1307674368000.0 },
66 { 16, 20922789888000.0 },
67 { 17, 355687428096000.0 },
68 { 18, 6402373705728000.0 },
69 { 19, 121645100408832000.0 },
70 { 20, 2432902008176640000.0 },
71 { 21, 51090942171709440000.0 },
72 { 22, 1124000727777607680000.0 },
73 { 23, 25852016738884976640000.0 },
74 { 24, 620448401733239439360000.0 },
75 { 25, 15511210043330985984000000.0 },
76 { 26, 403291461126605635584000000.0 },
77 { 27, 10888869450418352160768000000.0 },
78 { 28, 304888344611713860501504000000.0 },
79 { 29, 8841761993739701954543616000000.0 },
80 { 30, 265252859812191058636308480000000.0 },
81 { 31, 8222838654177922817725562880000000.0 },
82 { 32, 263130836933693530167218012160000000.0 },
83 { 33, 8683317618811886495518194401280000000.0 },
84 { 34, 2.95232799039604140847618609644e38 },
85 { 35, 1.03331479663861449296666513375e40 },
86 { 36, 3.71993326789901217467999448151e41 },
87 { 37, 1.37637530912263450463159795816e43 },
88 { 38, 5.23022617466601111760007224100e44 },
89 { 39, 2.03978820811974433586402817399e46 },
90 { 40, 8.15915283247897734345611269600e47 },
91 { 41, 3.34525266131638071081700620534e49 },
92 { 42, 1.40500611775287989854314260624e51 },
93 { 43, 6.04152630633738356373551320685e52 },
94 { 44, 2.65827157478844876804362581101e54 },
95 { 45, 1.19622220865480194561963161496e56 },
96 { 46, 5.50262215981208894985030542880e57 },
97 { 47, 2.58623241511168180642964355154e59 },
98 { 48, 1.24139155925360726708622890474e61 },
99 { 49, 6.08281864034267560872252163321e62 },
100 { 50, 3.04140932017133780436126081661e64 },
101 { 51, 1.55111875328738228022424301647e66 },
102 { 52, 8.06581751709438785716606368564e67 },
103 { 53, 4.27488328406002556429801375339e69 },
104 { 54, 2.30843697339241380472092742683e71 },
105 { 55, 1.26964033536582759259651008476e73 },
106 { 56, 7.10998587804863451854045647464e74 },
107 { 57, 4.05269195048772167556806019054e76 },
108 { 58, 2.35056133128287857182947491052e78 },
109 { 59, 1.38683118545689835737939019720e80 },
110 { 60, 8.32098711274139014427634118320e81 },
111 { 61, 5.07580213877224798800856812177e83 },
112 { 62, 3.14699732603879375256531223550e85 },
113 { 63, 1.982608315404440064116146708360e87 },
114 { 64, 1.268869321858841641034333893350e89 },
115 { 65, 8.247650592082470666723170306800e90 },
116 { 66, 5.443449390774430640037292402480e92 },
117 { 67, 3.647111091818868528824985909660e94 },
118 { 68, 2.480035542436830599600990418570e96 },
119 { 69, 1.711224524281413113724683388810e98 },
120 { 70, 1.197857166996989179607278372170e100 },
121 { 71, 8.504785885678623175211676442400e101 },
122 { 72, 6.123445837688608686152407038530e103 },
123 { 73, 4.470115461512684340891257138130e105 },
124 { 74, 3.307885441519386412259530282210e107 },
125 { 75, 2.480914081139539809194647711660e109 },
126 { 76, 1.885494701666050254987932260860e111 },
127 { 77, 1.451830920282858696340707840860e113 },
128 { 78, 1.132428117820629783145752115870e115 },
129 { 79, 8.946182130782975286851441715400e116 },
130 { 80, 7.156945704626380229481153372320e118 },
131 { 81, 5.797126020747367985879734231580e120 },
132 { 82, 4.753643337012841748421382069890e122 },
133 { 83, 3.945523969720658651189747118010e124 },
134 { 84, 3.314240134565353266999387579130e126 },
135 { 85, 2.817104114380550276949479442260e128 },
136 { 86, 2.422709538367273238176552320340e130 },
137 { 87, 2.107757298379527717213600518700e132 },
138 { 88, 1.854826422573984391147968456460e134 },
139 { 89, 1.650795516090846108121691926250e136 },
140 { 90, 1.485715964481761497309522733620e138 },
141 { 91, 1.352001527678402962551665687590e140 },
142 { 92, 1.243841405464130725547532432590e142 },
143 { 93, 1.156772507081641574759205162310e144 },
144 { 94, 1.087366156656743080273652852570e146 },
145 { 95, 1.032997848823905926259970209940e148 },
146 { 96, 9.916779348709496892095714015400e149 },
147 { 97, 9.619275968248211985332842594960e151 },
148 { 98, 9.426890448883247745626185743100e153 },
149 { 99, 9.332621544394415268169923885600e155 },
150 { 100,9.33262154439441526816992388563e157 },
151 { 101,9.42594775983835942085162312450e159 },
152 { 102,9.61446671503512660926865558700e161 },
153 { 103,9.90290071648618040754671525458e163 },
154 { 104,1.02990167451456276238485838648e166 },
155 { 105,1.08139675824029090050410130580e168 },
156 { 106,1.146280563734708354534347384148e170 },
157 { 107,1.226520203196137939351751701040e172 },
158 { 108,1.324641819451828974499891837120e174 },
159 { 109,1.443859583202493582204882102460e176 },
160 { 110,1.588245541522742940425370312710e178 },
161 { 111,1.762952551090244663872161047110e180 },
162 { 112,1.974506857221074023536820372760e182 },
163 { 113,2.231192748659813646596607021220e184 },
164 { 114,2.543559733472187557120132004190e186 },
165 { 115,2.925093693493015690688151804820e188 },
166 { 116,3.393108684451898201198256093590e190 },
167 { 117,3.96993716080872089540195962950e192 },
168 { 118,4.68452584975429065657431236281e194 },
169 { 119,5.57458576120760588132343171174e196 },
170 { 120,6.68950291344912705758811805409e198 },
171 { 121,8.09429852527344373968162284545e200 },
172 { 122,9.87504420083360136241157987140e202 },
173 { 123,1.21463043670253296757662432419e205 },
174 { 124,1.50614174151114087979501416199e207 },
175 { 125,1.88267717688892609974376770249e209 },
176 { 126,2.37217324288004688567714730514e211 },
177 { 127,3.01266001845765954480997707753e213 },
178 { 128,3.85620482362580421735677065923e215 },
179 { 129,4.97450422247728744039023415041e217 },
180 { 130,6.46685548922047367250730439554e219 },
181 { 131,8.47158069087882051098456875820e221 },
182 { 132,1.11824865119600430744996307608e224 },
183 { 133,1.48727070609068572890845089118e226 },
184 { 134,1.99294274616151887673732419418e228 },
185 { 135,2.69047270731805048359538766215e230 },
186 { 136,3.65904288195254865768972722052e232 },
187 { 137,5.01288874827499166103492629211e234 },
188 { 138,6.91778647261948849222819828311e236 },
189 { 139,9.61572319694108900419719561353e238 },
190 { 140,1.34620124757175246058760738589e241 },
191 { 141,1.89814375907617096942852641411e243 },
192 { 142,2.69536413788816277658850750804e245 },
193 { 143,3.85437071718007277052156573649e247 },
194 { 144,5.55029383273930478955105466055e249 },
195 { 145,8.04792605747199194484902925780e251 },
196 { 146,1.17499720439091082394795827164e254 },
197 { 147,1.72724589045463891120349865931e256 },
198 { 148,2.55632391787286558858117801578e258 },
199 { 149,3.80892263763056972698595524351e260 },
200 { 150,5.71338395644585459047893286526e262 },
201 { 151,8.62720977423324043162318862650e264 },
202 { 152,1.31133588568345254560672467123e267 },
203 { 153,2.00634390509568239477828874699e269 },
204 { 154,3.08976961384735088795856467036e271 },
205 { 155,4.78914290146339387633577523906e273 },
206 { 156,7.47106292628289444708380937294e275 },
207 { 157,1.17295687942641442819215807155e278 },
208 { 158,1.85327186949373479654360975305e280 },
209 { 159,2.94670227249503832650433950735e282 },
210 { 160,4.71472363599206132240694321176e284 },
211 { 161,7.59070505394721872907517857094e286 },
212 { 162,1.22969421873944943411017892849e289 },
213 { 163,2.00440157654530257759959165344e291 },
214 { 164,3.28721858553429622726333031164e293 },
215 { 165,5.42391066613158877498449501421e295 },
216 { 166,9.00369170577843736647426172359e297 },
217 { 167,1.50361651486499904020120170784e300 },
218 { 168,2.52607574497319838753801886917e302 },
219 { 169,4.26906800900470527493925188890e304 },
220 { 170,7.25741561530799896739672821113e306 },
221 { 171,1.24101807021766782342484052410e309L },
222 { 172,2.13455108077438865629072570146e311L },
223 { 173,3.69277336973969237538295546352e313L },
224 { 174,6.42542566334706473316634250653e315L },
225 { 175,1.12444949108573632830410993864e318L },
226 { 176,1.97903110431089593781523349201e320L },
227 { 177,3.50288505463028580993296328086e322L },
228 { 178,6.23513539724190874168067463993e324L },
229 { 179,1.11608923610630166476084076055e327L },
230 { 180,2.00896062499134299656951336898e329L },
231 { 181,3.63621873123433082379081919786e331L },
232 { 182,6.61791809084648209929929094011e333L },
233 { 183,1.21107901062490622417177024204e336L },
234 { 184,2.22838537954982745247605724535e338L },
235 { 185,4.12251295216718078708070590390e340L },
236 { 186,7.66787409103095626397011298130e342L },
237 { 187,1.43389245502278882136241112750e345L },
238 { 188,2.69571781544284298416133291969e347L },
239 { 189,5.09490667118697324006491921822e349L },
240 { 190,9.68032267525524915612334651460e351L },
241 { 191,1.84894163097375258881955918429e354L },
242 { 192,3.54996793146960497053355363384e356L },
243 { 193,6.85143810773633759312975851330e358L },
244 { 194,1.32917899290084949306717315158e361L },
245 { 195,2.59189903615665651148098764559e363L },
246 { 196,5.08012211086704676250273578535e365L },
247 { 197,1.00078405584080821221303894971e368L },
248 { 198,1.98155243056480026018181712043e370L },
249 { 199,3.94328933682395251776181606966e372L },
250 { 200,7.88657867364790503552363213932e374L }
251 };
252 /******************************************************************************/
253 /******************************************************************************/
GeocTipoCalcProd(void)254 int GeocTipoCalcProd(void)
255 {
256 //distingimos los tipos de cálculo
257 #if defined(CALCULO_PRODUCTO_MULT)
258 return GEOC_PROD_MULT;
259 #elif defined(CALCULO_PRODUCTO_LOG)
260 return GEOC_PROD_LOG;
261 #else
262 #error *****No se ha definido el método de cálculo de la función Producto()
263 #endif
264 }
265 /******************************************************************************/
266 /******************************************************************************/
Media(const double * datos,const size_t nDatos,const int inc)267 double Media(const double* datos,
268 const size_t nDatos,
269 const int inc)
270 {
271 //índice para recorrer un bucle
272 size_t i=0;
273 //posición del elemento de trabajo
274 size_t pos=0;
275 //variable para almacenar el resultado
276 double resultado=0.0;
277 ////////////////////////////////////////////////////////////////////////////
278 ////////////////////////////////////////////////////////////////////////////
279 //calculamos la posición del elemento inicial del vector
280 pos = GEOC_INICIO_VEC(nDatos,inc);
281 //recorremos el vector
282 for(i=0;i<nDatos;i++)
283 {
284 //vamos sumando los valores del vector
285 resultado += datos[pos];
286 //aumentamos el contador de posiciones
287 pos = inc>0 ? pos+(size_t)inc : pos-(size_t)abs(inc);
288 }
289 //calculamos la media
290 resultado /= (double)nDatos;
291 ////////////////////////////////////////////////////////////////////////////
292 ////////////////////////////////////////////////////////////////////////////
293 //salimos de la función
294 return resultado;
295 }
296 /******************************************************************************/
297 /******************************************************************************/
Varianza(const double * datos,const size_t nDatos,const int inc,const double media)298 double Varianza(const double* datos,
299 const size_t nDatos,
300 const int inc,
301 const double media)
302 {
303 //índice para recorrer un bucle
304 size_t i=0;
305 //posición del elemento de trabajo
306 size_t pos=0;
307 //residuo
308 double res=0.0;
309 //variable para almacenar el resultado
310 double resultado=0.0;
311 ////////////////////////////////////////////////////////////////////////////
312 ////////////////////////////////////////////////////////////////////////////
313 //calculamos la posición del elemento inicial del vector
314 pos = GEOC_INICIO_VEC(nDatos,inc);
315 //recorremos el vector
316 for(i=0;i<nDatos;i++)
317 {
318 //calculamos el residuo
319 res = datos[pos]-media;
320 //vamos sumando los cuadrados del valor menos la media
321 resultado += res*res;
322 //aumentamos el contador de posiciones
323 pos = inc>0 ? pos+(size_t)inc : pos-(size_t)abs(inc);
324 }
325 //calculamos la varianza
326 resultado /= (double)(nDatos-1);
327 ////////////////////////////////////////////////////////////////////////////
328 ////////////////////////////////////////////////////////////////////////////
329 //salimos de la función
330 return resultado;
331 }
332 /******************************************************************************/
333 /******************************************************************************/
MediaPonderada(const double * datos,const size_t nDatos,const int incDatos,const double * pesos,const int incPesos)334 double MediaPonderada(const double* datos,
335 const size_t nDatos,
336 const int incDatos,
337 const double* pesos,
338 const int incPesos)
339 {
340 //índice para recorrer un bucle
341 size_t i=0;
342 //posiciones de los elementos de trabajo
343 size_t posDatos=0,posPesos;
344 //variable auxiliar
345 double sumaPesos=0.0;
346 //variable para almacenar el resultado
347 double resultado=0.0;
348 ////////////////////////////////////////////////////////////////////////////
349 ////////////////////////////////////////////////////////////////////////////
350 //calculamos la posición del elemento inicial de los vectores
351 posDatos = GEOC_INICIO_VEC(nDatos,incDatos);
352 posPesos = GEOC_INICIO_VEC(nDatos,incPesos);
353 //recorremos el vector
354 for(i=0;i<nDatos;i++)
355 {
356 //vamos sumando los valores del vector
357 resultado += datos[posDatos]*pesos[posPesos];
358 //aumentamos el contador de posiciones
359 posDatos = incDatos>0 ? posDatos+(size_t)incDatos
360 : posDatos-(size_t)abs(incDatos);
361 posPesos = incPesos>0 ? posPesos+(size_t)incPesos
362 : posPesos-(size_t)abs(incPesos);
363 //vamos sumando los pesos
364 sumaPesos += pesos[posPesos];
365 }
366 //calculamos la media
367 resultado /= sumaPesos;
368 ////////////////////////////////////////////////////////////////////////////
369 ////////////////////////////////////////////////////////////////////////////
370 //salimos de la función
371 return resultado;
372 }
373 /******************************************************************************/
374 /******************************************************************************/
Mediana(const double * datos,const size_t nDatos,const int inc)375 double Mediana(const double* datos,
376 const size_t nDatos,
377 const int inc)
378 {
379 //posición del elemento central
380 size_t pos=0;
381 //posición de inicio del vector
382 size_t posInicio;
383 //variable auxiliar
384 size_t posAux=0;
385 //variable para almacenar el resultado
386 double resultado=0.0;
387 ////////////////////////////////////////////////////////////////////////////
388 ////////////////////////////////////////////////////////////////////////////
389 //calculamos la posición del elemento inicial del vector
390 posInicio = GEOC_INICIO_VEC(nDatos,inc);
391 //distinguimos entre número par o impar de datos
392 if(nDatos%2)
393 {
394 //el número de datos es impar, calculamos la posición del punto medio
395 pos = (size_t)ceil((double)nDatos/2.0)-1;
396 //calculamos la posición en memoria
397 pos = inc>0 ? posInicio+pos*(size_t)inc
398 : posInicio-pos*(size_t)abs(inc);
399 //extraemos el dato
400 resultado = datos[pos];
401 }
402 else
403 {
404 //el número de datos es par, calculamos la posición de los puntos medios
405 pos = (size_t)((double)nDatos/2.0)-1;
406 //calculamos la posición en memoria
407 pos = inc>0 ? posInicio+pos*(size_t)inc
408 : posInicio-pos*(size_t)abs(inc);
409 //calculamos la posición del siguiente elemento del vector
410 posAux = inc>0 ? pos+(size_t)inc : pos-(size_t)abs(inc);
411 //calculamos la mediana
412 resultado = (datos[pos]+datos[posAux])/2.0;
413 }
414 ////////////////////////////////////////////////////////////////////////////
415 ////////////////////////////////////////////////////////////////////////////
416 //salimos de la función
417 return resultado;
418 }
419 /******************************************************************************/
420 /******************************************************************************/
ProductoMult(size_t inicio,size_t fin)421 double ProductoMult(size_t inicio,
422 size_t fin)
423 {
424 //variables intermedias
425 size_t menor=0.0,mayor=0.0;
426 //variable para almacenar el resultado
427 double resultado=1.0;
428 ////////////////////////////////////////////////////////////////////////////
429 ////////////////////////////////////////////////////////////////////////////
430 //resolvemos el problema de forma iterativa
431 if((inicio==0)||(fin==0))
432 {
433 resultado = 0.0;
434 }
435 else
436 {
437 //seleccionamos los elementos mayor y menor
438 menor = inicio<=fin ? inicio : fin;
439 mayor = inicio>fin ? inicio : fin;
440 //multiplicamos
441 while(menor<=mayor)
442 {
443 resultado *= (double)menor;
444 menor++;
445 }
446 }
447 ////////////////////////////////////////////////////////////////////////////
448 ////////////////////////////////////////////////////////////////////////////
449 //salimos de la función
450 return resultado;
451 }
452 /******************************************************************************/
453 /******************************************************************************/
ProductoMultLD(size_t inicio,size_t fin)454 long double ProductoMultLD(size_t inicio,
455 size_t fin)
456 {
457 //variables intermedias
458 size_t menor=0.0,mayor=0.0;
459 //variable para almacenar el resultado
460 long double resultado=1.0;
461 ////////////////////////////////////////////////////////////////////////////
462 ////////////////////////////////////////////////////////////////////////////
463 //resolvemos el problema de forma iterativa
464 if((inicio==0)||(fin==0))
465 {
466 resultado = 0.0;
467 }
468 else
469 {
470 //seleccionamos los elementos mayor y menor
471 menor = inicio<=fin ? inicio : fin;
472 mayor = inicio>fin ? inicio : fin;
473 //multiplicamos
474 while(menor<=mayor)
475 {
476 resultado *= (long double)menor;
477 menor++;
478 }
479 }
480 ////////////////////////////////////////////////////////////////////////////
481 ////////////////////////////////////////////////////////////////////////////
482 //salimos de la función
483 return resultado;
484 }
485 /******************************************************************************/
486 /******************************************************************************/
ProductoLog(size_t inicio,size_t fin)487 double ProductoLog(size_t inicio,
488 size_t fin)
489 {
490 //variables intermedias
491 size_t menor=0.0,mayor=0.0;
492 //variable para almacenar el resultado
493 double resultado=0.0;
494 ////////////////////////////////////////////////////////////////////////////
495 ////////////////////////////////////////////////////////////////////////////
496 //resolvemos el problema de forma iterativa
497 if((inicio==0)||(fin==0))
498 {
499 resultado = 0.0;
500 }
501 else
502 {
503 //seleccionamos los elementos mayor y menor
504 menor = inicio<=fin ? inicio : fin;
505 mayor = inicio>fin ? inicio : fin;
506 //multiplicamos
507 while(menor<=mayor)
508 {
509 resultado += log((double)menor);
510 menor++;
511 }
512 }
513 ////////////////////////////////////////////////////////////////////////////
514 ////////////////////////////////////////////////////////////////////////////
515 //salimos de la función
516 return round(pow(GEOC_CONST_E,resultado));
517 }
518 /******************************************************************************/
519 /******************************************************************************/
ProductoLogLD(size_t inicio,size_t fin)520 long double ProductoLogLD(size_t inicio,
521 size_t fin)
522 {
523 //variables intermedias
524 size_t menor=0.0,mayor=0.0;
525 //variable para almacenar el resultado
526 long double resultado=0.0;
527 ////////////////////////////////////////////////////////////////////////////
528 ////////////////////////////////////////////////////////////////////////////
529 //resolvemos el problema de forma iterativa
530 if((inicio==0)||(fin==0))
531 {
532 resultado = 0.0;
533 }
534 else
535 {
536 //seleccionamos los elementos mayor y menor
537 menor = inicio<=fin ? inicio : fin;
538 mayor = inicio>fin ? inicio : fin;
539 //multiplicamos
540 while(menor<=mayor)
541 {
542 resultado += (long double)log((double)menor);
543 menor++;
544 }
545 }
546 ////////////////////////////////////////////////////////////////////////////
547 ////////////////////////////////////////////////////////////////////////////
548 //salimos de la función
549 return roundl(powl((long double)GEOC_CONST_E,resultado));
550 }
551 /******************************************************************************/
552 /******************************************************************************/
Producto(size_t inicio,size_t fin)553 double Producto(size_t inicio,
554 size_t fin)
555 {
556 //variable para almacenar el resultado
557 double resultado=0;
558 ////////////////////////////////////////////////////////////////////////////
559 ////////////////////////////////////////////////////////////////////////////
560 //método de cálculo para números mayores que GEOC_MATE_CONST_NMAXFAC
561 #if defined(CALCULO_PRODUCTO_MULT)
562 resultado = ProductoMult(inicio,fin);
563 #elif defined(CALCULO_PRODUCTO_LOG)
564 resultado = ProductoLog(inicio,fin);
565 #else
566 #error *****No se ha definido el método de cálculo de la función Producto()
567 #endif
568 ////////////////////////////////////////////////////////////////////////////
569 ////////////////////////////////////////////////////////////////////////////
570 //salimos de la función
571 return resultado;
572 }
573 /******************************************************************************/
574 /******************************************************************************/
ProductoLD(size_t inicio,size_t fin)575 long double ProductoLD(size_t inicio,
576 size_t fin)
577 {
578 //variable para almacenar el resultado
579 long double resultado=0;
580 ////////////////////////////////////////////////////////////////////////////
581 ////////////////////////////////////////////////////////////////////////////
582 //método de cálculo para números mayores que GEOC_MATE_CONST_NMAXFAC
583 #if defined(CALCULO_PRODUCTO_MULT)
584 resultado = ProductoMultLD(inicio,fin);
585 #elif defined(CALCULO_PRODUCTO_LOG)
586 resultado = ProductoLogLD(inicio,fin);
587 #else
588 #error *****No se ha definido el método de cálculo de la función Producto()
589 #endif
590 ////////////////////////////////////////////////////////////////////////////
591 ////////////////////////////////////////////////////////////////////////////
592 //salimos de la función
593 return resultado;
594 }
595 /******************************************************************************/
596 /******************************************************************************/
FactorialMult(size_t numero)597 double FactorialMult(size_t numero)
598 {
599 //variable para almacenar el resultado
600 double resultado=1.0;
601 ////////////////////////////////////////////////////////////////////////////
602 ////////////////////////////////////////////////////////////////////////////
603 //resolvemos el problema
604 if(numero)
605 {
606 resultado = ProductoMult((size_t)1,numero);
607 }
608 ////////////////////////////////////////////////////////////////////////////
609 ////////////////////////////////////////////////////////////////////////////
610 //salimos de la función
611 return resultado;
612 }
613 /******************************************************************************/
614 /******************************************************************************/
FactorialMultLD(size_t numero)615 long double FactorialMultLD(size_t numero)
616 {
617 //variable para almacenar el resultado
618 long double resultado=1.0;
619 ////////////////////////////////////////////////////////////////////////////
620 ////////////////////////////////////////////////////////////////////////////
621 //resolvemos el problema
622 if(numero)
623 {
624 resultado = ProductoMultLD((size_t)1,numero);
625 }
626 ////////////////////////////////////////////////////////////////////////////
627 ////////////////////////////////////////////////////////////////////////////
628 //salimos de la función
629 return resultado;
630 }
631 /******************************************************************************/
632 /******************************************************************************/
FactorialLog(size_t numero)633 double FactorialLog(size_t numero)
634 {
635 //variable para almacenar el resultado
636 double resultado=1.0;
637 ////////////////////////////////////////////////////////////////////////////
638 ////////////////////////////////////////////////////////////////////////////
639 //resolvemos el problema
640 if(numero)
641 {
642 resultado = ProductoLog((size_t)1,numero);
643 }
644 ////////////////////////////////////////////////////////////////////////////
645 ////////////////////////////////////////////////////////////////////////////
646 //salimos de la función
647 return resultado;
648 }
649 /******************************************************************************/
650 /******************************************************************************/
FactorialLogLD(size_t numero)651 long double FactorialLogLD(size_t numero)
652 {
653 //variable para almacenar el resultado
654 long double resultado=1.0;
655 ////////////////////////////////////////////////////////////////////////////
656 ////////////////////////////////////////////////////////////////////////////
657 //resolvemos el problema
658 if(numero)
659 {
660 resultado = ProductoLogLD((size_t)1,numero);
661 }
662 ////////////////////////////////////////////////////////////////////////////
663 ////////////////////////////////////////////////////////////////////////////
664 //salimos de la función
665 return resultado;
666 }
667 /******************************************************************************/
668 /******************************************************************************/
Factorial(size_t numero)669 double Factorial(size_t numero)
670 {
671 //variable para almacenar el resultado
672 double resultado=0;
673 ////////////////////////////////////////////////////////////////////////////
674 ////////////////////////////////////////////////////////////////////////////
675 //calculamos el factorial
676 if(numero<=GEOC_MATE_CONST_DBL_NMAXFAC)
677 {
678 resultado = (double)__tablaFacExpl[numero].valor;
679 }
680 else
681 {
682 //método de cálculo para números mayores que GEOC_MATE_CONST_DBL_NMAXFAC
683 #if defined(CALCULO_PRODUCTO_MULT)
684 resultado = FactorialMult(numero);
685 #elif defined(CALCULO_PRODUCTO_LOG)
686 resultado = FactorialLog(numero);
687 #else
688 #error *****No se ha definido el método de cálculo de la función Factorial()
689 #endif
690 }
691 ////////////////////////////////////////////////////////////////////////////
692 ////////////////////////////////////////////////////////////////////////////
693 //salimos de la función
694 return resultado;
695 }
696 /******************************************************************************/
697 /******************************************************************************/
FactorialLD(size_t numero)698 long double FactorialLD(size_t numero)
699 {
700 //variable para almacenar el resultado
701 long double resultado=0;
702 ////////////////////////////////////////////////////////////////////////////
703 ////////////////////////////////////////////////////////////////////////////
704 //calculamos el factorial
705 if(numero<=GEOC_MATE_CONST_LDBL_NMAXFAC)
706 {
707 resultado = __tablaFacExpl[numero].valor;
708 }
709 else
710 {
711 //cálculo para números mayores que GEOC_MATE_CONST_LDBL_NMAXFAC
712 #if defined(CALCULO_PRODUCTO_MULT)
713 resultado = FactorialMultLD(numero);
714 #elif defined(CALCULO_PRODUCTO_LOG)
715 resultado = FactorialLogLD(numero);
716 #else
717 #error *****No se ha definido el método de cálculo de la función Factorial()
718 #endif
719 }
720 ////////////////////////////////////////////////////////////////////////////
721 ////////////////////////////////////////////////////////////////////////////
722 //salimos de la función
723 return resultado;
724 }
725 /******************************************************************************/
726 /******************************************************************************/
ProductoVectorial(const double x1,const double y1,const double z1,const double x2,const double y2,const double z2,double * x,double * y,double * z)727 void ProductoVectorial(const double x1,
728 const double y1,
729 const double z1,
730 const double x2,
731 const double y2,
732 const double z2,
733 double* x,
734 double* y,
735 double* z)
736 {
737 //calculamos las componentes
738 *x = y1*z2-y2*z1;
739 *y = x2*z1-x1*z2;
740 *z = x1*y2-x2*y1;
741 ////////////////////////////////////////////////////////////////////////////
742 ////////////////////////////////////////////////////////////////////////////
743 //salimos de la función
744 return;
745 }
746 /******************************************************************************/
747 /******************************************************************************/
SinCosRecurrencia(const double anguloIni,const double incAngulo,const size_t numValores,double * seno,const size_t incSeno,double * coseno,const size_t incCoseno)748 void SinCosRecurrencia(const double anguloIni,
749 const double incAngulo,
750 const size_t numValores,
751 double* seno,
752 const size_t incSeno,
753 double* coseno,
754 const size_t incCoseno)
755 {
756 //índice para recorrer bucles
757 size_t i=0;
758 //posiciones en los vectores de trabajo
759 size_t pS=0,pC=0;
760 //variables auxiliares
761 double aux=sin(incAngulo/2.0);
762 double alfa=2.0*aux*aux,beta=sin(incAngulo);
763 ////////////////////////////////////////////////////////////////////////////
764 ////////////////////////////////////////////////////////////////////////////
765 //seno y coseno del primero elemento de la serie
766 seno[pS] = sin(anguloIni);
767 coseno[pC] = cos(anguloIni);
768 //recorremos el resto de valores
769 for(i=1;i<numValores;i++)
770 {
771 //calculamos las posiciones en los vectores de salida
772 pS += incSeno;
773 pC += incCoseno;
774 //vamos calculando los senos y los cosenos
775 seno[pS] = seno[pS-incSeno]-
776 (alfa*seno[pS-incSeno]-beta*coseno[pC-incCoseno]);
777 coseno[pC] = coseno[pC-incCoseno]-
778 (alfa*coseno[pC-incCoseno]+beta*seno[pS-incSeno]);
779 }
780 ////////////////////////////////////////////////////////////////////////////
781 ////////////////////////////////////////////////////////////////////////////
782 //salimos de la función
783 return;
784 }
785 /******************************************************************************/
786 /******************************************************************************/
SplineCubicoNatural(double * x,double * y,size_t nDatos,double * xInterp,size_t nInterp)787 void SplineCubicoNatural(double* x,
788 double* y,
789 size_t nDatos,
790 double* xInterp,
791 size_t nInterp)
792 {
793 //sentencias para que no salgan avisos en la compilación
794 x = x;
795 y = y;
796 nDatos = nDatos;
797 xInterp = xInterp;
798 nInterp = nInterp;
799 ////////////////////////////////////////////////////////////////////////////
800 ////////////////////////////////////////////////////////////////////////////
801 //salimos de la función
802 return;
803 }
804 /******************************************************************************/
805 /******************************************************************************/
806 /** @} */
807 /******************************************************************************/
808 /******************************************************************************/
809 /* kate: encoding utf-8; end-of-line unix; syntax c; indent-mode cstyle; */
810 /* kate: replace-tabs on; space-indent on; tab-indents off; indent-width 4; */
811 /* kate: line-numbers on; folding-markers on; remove-trailing-space on; */
812 /* kate: backspace-indents on; show-tabs on; */
813 /* kate: word-wrap-column 80; word-wrap-marker-color #D2D2D2; word-wrap off; */
814