Bug Summary

File:OMCompiler/SimulationRuntime/c/./simulation/solver/kinsolSolver.c
Warning:line 349, column 17
Value stored to 'threadData' during its initialization is never read

Annotated Source Code

[?] Use j/k keys for keyboard navigation

1/*
2 * This file is part of OpenModelica.
3 *
4 * Copyright (c) 1998-2014, Open Source Modelica Consortium (OSMC),
5 * c/o Linköpings universitet, Department of Computer and Information Science,
6 * SE-58183 Linköping, Sweden.
7 *
8 * All rights reserved.
9 *
10 * THIS PROGRAM IS PROVIDED UNDER THE TERMS OF THE BSD NEW LICENSE OR THE
11 * GPL VERSION 3 LICENSE OR THE OSMC PUBLIC LICENSE (OSMC-PL) VERSION 1.2.
12 * ANY USE, REPRODUCTION OR DISTRIBUTION OF THIS PROGRAM CONSTITUTES
13 * RECIPIENT'S ACCEPTANCE OF THE OSMC PUBLIC LICENSE OR THE GPL VERSION 3,
14 * ACCORDING TO RECIPIENTS CHOICE.
15 *
16 * The OpenModelica software and the OSMC (Open Source Modelica Consortium)
17 * Public License (OSMC-PL) are obtained from OSMC, either from the above
18 * address, from the URLs: http://www.openmodelica.org or
19 * http://www.ida.liu.se/projects/OpenModelica, and in the OpenModelica
20 * distribution. GNU version 3 is obtained from:
21 * http://www.gnu.org/copyleft/gpl.html. The New BSD License is obtained from:
22 * http://www.opensource.org/licenses/BSD-3-Clause.
23 *
24 * This program is distributed WITHOUT ANY WARRANTY; without even the implied
25 * warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE, EXCEPT AS
26 * EXPRESSLY SET FORTH IN THE BY RECIPIENT SELECTED SUBSIDIARY LICENSE
27 * CONDITIONS OF OSMC-PL.
28 *
29 */
30
31/*! \file kinsolSolver.c
32 */
33
34#include "omc_config.h"
35#include "nonlinearSystem.h"
36#include "kinsolSolver.h"
37#include "simulation/simulation_info_json.h"
38#include "simulation/options.h"
39#include "util/omc_error.h"
40#include "omc_math.h"
41
42#ifdef WITH_SUNDIALS
43
44/* adrpo: on mingw link with static sundials */
45#if defined(__MINGW32__)
46#define LINK_SUNDIALS_STATIC
47#endif
48
49#include <math.h>
50#include <stdlib.h>
51#include <string.h> /* memcpy */
52
53#include "util/varinfo.h"
54#include "openmodelica.h"
55#include "openmodelica_func.h"
56#include "model_help.h"
57#include "util/read_matlab4.h"
58#include "events.h"
59
60#include <string.h>
61#include <stdio.h>
62#include <stdlib.h>
63#include <math.h>
64
65#include <kinsol/kinsol.h>
66#include <kinsol/kinsol_impl.h>
67#include <kinsol/kinsol_dense.h>
68#include <kinsol/kinsol_lapack.h>
69#include <kinsol/kinsol_direct.h>
70#include <kinsol/kinsol_klu.h>
71#include <sundials/sundials_nvector.h>
72#include <nvector/nvector_serial.h>
73#include <sundials/sundials_types.h>
74#include <sundials/sundials_math.h>
75
76
77/* constants */
78#define RETRY_MAX5 5
79#define FTOL_WITH_LESS_ACCURANCY1.e-6 1.e-6
80
81/* readability */
82#define INITIAL_EXTRAPOLATION0 0
83#define INITIAL_OLDVALUES1 1
84
85#define SCALING_NOMINALSTART1 1
86#define SCALING_ONES2 2
87#define SCALING_JACOBIAN3 3
88
89
90
91typedef struct NLS_KINSOL_USERDATA
92{
93 DATA* data;
94 threadData_t* threadData;
95
96 int sysNumber;
97}NLS_KINSOL_USERDATA;
98
99typedef struct NLS_KINSOL_DATA
100{
101 /* ### configuration ### */
102 int linearSolverMethod; /* specifies the method to solve the underlying linear problem */
103 int nonLinearSystemNumber;
104 int kinsolStrategy;
105 int retries;
106 int solved; /* if the system is once solved reuse linear matrix information */
107 int nominalJac;
108
109 /* ### tolerances ### */
110 double fnormtol; /* function tolerance */
111 double scsteptol; /* step tolerance */
112 double maxstepfactor; /* maximum newton step factor mxnewtstep = maxstepfactor * norm2(xScaling) */
113
114 /* ### work arrays ### */
115 N_Vector initialGuess;
116 N_Vector xScale;
117 N_Vector fScale;
118 N_Vector fRes;
119 N_Vector fTmp;
120
121 int iflag;
122 long countResCalls; /* case of sparse function not avaiable */
123
124 /* ### kinsol internal data */
125 void* kinsolMemory;
126 NLS_KINSOL_USERDATA userData;
127
128 /* settings */
129 int size;
130 int nnz;
131
132}NLS_KINSOL_DATA;
133
134static int nlsKinsolResiduals(N_Vector z, N_Vector f, void *userData);
135static void nlsKinsolErrorPrint(int error_code, const char *module, const char *function, char *msg, void *userData);
136static void nlsKinsolInfoPrint(const char *module, const char *function, char *msg, void *userData);
137static int nlsSparseJac(N_Vector x, N_Vector fx, SlsMat Jac, void *userData, N_Vector tmp1, N_Vector tmp2);
138static int nlsSparseSymJac(N_Vector x, N_Vector fx, SlsMat Jac, void *userData, N_Vector tmp1, N_Vector tmp2);
139static int nlsDenseJac(long int N, N_Vector x, N_Vector fx, DlsMat Jac, void *userData, N_Vector tmp1, N_Vector tmp2);
140static void nlsKinsolJacSumSparse(SlsMat mat);
141static void nlsKinsolJacSumDense(DlsMat mat);
142
143
144
145int checkReturnFlag(int flag)
146{
147 int retVal = flag;
148 switch(flag)
149 {
150 case KIN_SUCCESS0:
151 retVal = 0;
152 break;
153 case KIN_MEM_NULL-1:
154 case KIN_MEM_FAIL-4:
155 case KIN_ILL_INPUT-2:
156 retVal = -1;
157 break;
158 }
159 return retVal;
160}
161
162static
163void nlsKinsolConfigSetup(NLS_KINSOL_DATA *kinsolData)
164{
165 /* configuration */
166 KINSetFuncNormTol(kinsolData->kinsolMemory, kinsolData->fnormtol);
167 KINSetScaledStepTol(kinsolData->kinsolMemory, kinsolData->scsteptol);
168 KINSetNumMaxIters(kinsolData->kinsolMemory, 100*kinsolData->size);
169 kinsolData->kinsolStrategy = KIN_LINESEARCH1;
170 /* configuration for exact Newton */
171 /*
172 KINSetMaxSetupCalls(kinsolData->kinsolMemory, 1);
173 KINSetMaxSubSetupCalls(kinsolData->kinsolMemory, 1);
174 */
175
176 KINSetNoInitSetup(kinsolData->kinsolMemory, FALSE0);
177
178 kinsolData->retries = 0;
179 kinsolData->countResCalls = 0;
180}
181
182static void resetKinsolMemory(NLS_KINSOL_DATA *kinsolData, NONLINEAR_SYSTEM_DATA *nlsData)
183{
184 int flag, size = kinsolData->size, printLevel;
185 if (kinsolData->kinsolMemory) {
186 KINFree((void*)&kinsolData->kinsolMemory);
187 }
188 kinsolData->kinsolMemory = KINCreate();
189 /* setup user defined functions */
190 KINSetErrHandlerFn(kinsolData->kinsolMemory, nlsKinsolErrorPrint, kinsolData);
191 KINSetInfoHandlerFn(kinsolData->kinsolMemory, nlsKinsolInfoPrint, kinsolData);
192 KINSetUserData(kinsolData->kinsolMemory, (void*)&(kinsolData->userData));
193 flag = KINInit(kinsolData->kinsolMemory, nlsKinsolResiduals, kinsolData->initialGuess);
194 if (checkReturnFlag(flag)){
195 errorStreamPrint(LOG_STDOUT, 0, "##KINSOL## Something goes wrong while initialize KINSOL solver!");
196 }
197 /* Specify linear solver and/or corresponding jacobian function*/
198 if (kinsolData->linearSolverMethod == NLS_LS_KLU)
199 {
200 if(nlsData->isPatternAvailable)
201 {
202 kinsolData->nnz = nlsData->sparsePattern->numberOfNoneZeros;
203 flag = KINKLU(kinsolData->kinsolMemory, size, kinsolData->nnz);
204 if (checkReturnFlag(flag)){
205 errorStreamPrint(LOG_STDOUT, 0, "##KINSOL## Something goes wrong while initialize KINSOL solver!");
206 }
207 if (nlsData->analyticalJacobianColumn != NULL((void*)0)){
208 flag = KINSlsSetSparseJacFn(kinsolData->kinsolMemory, nlsSparseSymJac);
209 }
210 else {
211 flag = KINSlsSetSparseJacFn(kinsolData->kinsolMemory, nlsSparseJac);
212 }
213 if (checkReturnFlag(flag)){
214 errorStreamPrint(LOG_STDOUT, 0, "##KINSOL## Something goes wrong while initialize KINSOL Sparse Solver!");
215 }
216 } else {
217 flag = KINDense(kinsolData->kinsolMemory, size);
218 if (checkReturnFlag(flag)){
219 errorStreamPrint(LOG_STDOUT, 0, "##KINSOL## Something goes wrong while initialize KINSOL solver!");
220 }
221 }
222 }
223 else if (kinsolData->linearSolverMethod == NLS_LS_TOTALPIVOT)
224 {
225 flag = KINDense(kinsolData->kinsolMemory, size);
226 if (checkReturnFlag(flag)){
227 errorStreamPrint(LOG_STDOUT, 0, "##KINSOL## Something goes wrong while initialize KINSOL solver!");
228 }
229 }
230 else if (kinsolData->linearSolverMethod == NLS_LS_LAPACK)
231 {
232 flag = KINDense(kinsolData->kinsolMemory, size);
233 if (checkReturnFlag(flag)){
234 errorStreamPrint(LOG_STDOUT, 0, "##KINSOL## Something goes wrong while initialize KINSOL solver!");
235 }
236 flag = KINDlsSetDenseJacFn(kinsolData->kinsolMemory, nlsDenseJac);
237 if (checkReturnFlag(flag)){
238 errorStreamPrint(LOG_STDOUT, 0, "##KINSOL## Something goes wrong while initialize KINSOL Sparse Solver!");
239 }
240 }
241
242 /* configuration */
243 nlsKinsolConfigSetup(kinsolData);
244
245 /* debug print level of kinsol */
246 if (ACTIVE_STREAM(LOG_NLS_V)(useStream[LOG_NLS_V])) {
247 printLevel = 3;
248 } else if (ACTIVE_STREAM(LOG_NLS)(useStream[LOG_NLS])) {
249 printLevel = 1;
250 } else {
251 printLevel = 0;
252 }
253 KINSetPrintLevel(kinsolData->kinsolMemory, printLevel);
254}
255
256int nlsKinsolAllocate(int size, NONLINEAR_SYSTEM_DATA *nlsData, int linearSolverMethod)
257{
258 NLS_KINSOL_DATA *kinsolData = (NLS_KINSOL_DATA*) malloc(sizeof(NLS_KINSOL_DATA));
259
260 /* allocate system data */
261 nlsData->solverData = (void*)kinsolData;
262
263 kinsolData->size = size;
264 kinsolData->linearSolverMethod = linearSolverMethod;
265 kinsolData->solved = 0;
266
267 kinsolData->fnormtol = newtonFTol; /* function tolerance */
268 kinsolData->scsteptol = newtonXTol; /* step tolerance */
269
270 kinsolData->maxstepfactor = maxStepFactor; /* step tolerance */
271 kinsolData->nominalJac = 0; /* calculate for scaling the scaled matrix */
272
273 kinsolData->initialGuess = N_VNew_Serial(size);
274 kinsolData->xScale = N_VNew_Serial(size);
275 kinsolData->fScale = N_VNew_Serial(size);
276 kinsolData->fRes = N_VNew_Serial(size);
277 kinsolData->fTmp = N_VNew_Serial(size);
278
279 kinsolData->kinsolMemory = NULL((void*)0);
280
281 resetKinsolMemory(kinsolData, nlsData);
282
283 return 0;
284}
285
286int nlsKinsolFree(void** solverData)
287{
288 NLS_KINSOL_DATA* kinsolData = (NLS_KINSOL_DATA*) *solverData;
289
290 KINFree((void*)&kinsolData->kinsolMemory);
291
292 N_VDestroy_Serial(kinsolData->initialGuess);
293 N_VDestroy_Serial(kinsolData->xScale);
294 N_VDestroy_Serial(kinsolData->fScale);
295 N_VDestroy_Serial(kinsolData->fRes);
296 N_VDestroy_Serial(kinsolData->fTmp);
297 free(kinsolData);
298
299 return 0;
300}
301
302 /*! \fn nlsKinsolResiduals
303 *
304 * \param [in] [x]
305 * \param [out] [f]
306 * \param [ref] [user_data]
307 *
308 */
309static int nlsKinsolResiduals(N_Vector x, N_Vector f, void *userData)
310{
311 double *xdata = NV_DATA_S(x)( ( (N_VectorContent_Serial)(x->content) )->data );
312 double *fdata = NV_DATA_S(f)( ( (N_VectorContent_Serial)(f->content) )->data );
313
314 NLS_KINSOL_USERDATA *kinsolUserData = (NLS_KINSOL_USERDATA*) userData;
315 DATA* data = kinsolUserData->data;
316 threadData_t *threadData = kinsolUserData->threadData;
317 int sysNumber = kinsolUserData->sysNumber;
318 void *dataAndThreadData[2] = {data, threadData};
319 NONLINEAR_SYSTEM_DATA *nlsData = &(data->simulationInfo->nonlinearSystemData[sysNumber]);
320 NLS_KINSOL_DATA* kinsolData = (NLS_KINSOL_DATA*) nlsData->solverData;
321 long eqSystemNumber = nlsData->equationIndex;
322 int iflag = 1;
323
324 kinsolData->countResCalls++;
325
326#ifndef OMC_EMCC
327 MMC_TRY_INTERNAL(simulationJumpBuffer){ jmp_buf new_mmc_jumper, *old_jumper __attribute__((unused))
= threadData->simulationJumpBuffer; threadData->simulationJumpBuffer
= &new_mmc_jumper; if (_setjmp (new_mmc_jumper) == 0) {
328#endif
329
330 /* call residual function */
331 data->simulationInfo->nonlinearSystemData[sysNumber].residualFunc(dataAndThreadData, xdata, fdata, (const int*) &iflag);
332 iflag = 0;
333
334#ifndef OMC_EMCC
335 MMC_CATCH_INTERNAL(simulationJumpBuffer)} threadData->simulationJumpBuffer = old_jumper;mmc_catch_dummy_fn
();}
336#endif
337
338 return iflag;
339}
340
341/*
342 * function calculates a jacobian matrix
343 */
344static
345int nlsDenseJac(long int N, N_Vector vecX, N_Vector vecFX, DlsMat Jac, void *userData, N_Vector tmp1, N_Vector tmp2)
346{
347 NLS_KINSOL_USERDATA *kinsolUserData = (NLS_KINSOL_USERDATA*) userData;
348 DATA* data = kinsolUserData->data;
349 threadData_t *threadData = kinsolUserData->threadData;
Value stored to 'threadData' during its initialization is never read
350 int sysNumber = kinsolUserData->sysNumber;
351 NONLINEAR_SYSTEM_DATA *nlsData = &(data->simulationInfo->nonlinearSystemData[sysNumber]);
352 NLS_KINSOL_DATA* kinsolData = (NLS_KINSOL_DATA*) nlsData->solverData;
353
354 /* prepare variables */
355 double *x = N_VGetArrayPointer(vecX);
356 double *fx = N_VGetArrayPointer(vecFX);
357 double *xScaling = NV_DATA_S(kinsolData->xScale)( ( (N_VectorContent_Serial)(kinsolData->xScale->content
) )->data )
;
358 double *fRes = NV_DATA_S(kinsolData->fRes)( ( (N_VectorContent_Serial)(kinsolData->fRes->content)
)->data )
;
359 double xsave, xscale, sign;
360 double delta_hh;
361 const double delta_h = sqrt(DBL_EPSILON2.2204460492503131e-16*2e1);
362
363 long int i,j;
364
365 /* performance measurement */
366 rt_ext_tp_tick(&nlsData->jacobianTimeClock);
367
368 for(i = 0; i < N; i++)
369 {
370 xsave = x[i];
371 delta_hh = delta_h * (fabs(xsave) + 1.0);
372 if ((xsave + delta_hh >= nlsData->max[i]))
373 delta_hh *= -1;
374 x[i] += delta_hh;
375
376 /* Calculate difference quotient */
377 nlsKinsolResiduals(vecX, kinsolData->fRes, userData);
378
379 /* Calculate scaled difference quotient */
380 delta_hh = 1. / delta_hh;
381
382 for(j = 0; j < N; j++)
383 {
384 if (kinsolData->nominalJac){
385 DENSE_ELEM(Jac, j, i)((Jac->cols)[i][j]) = (fRes[j] - fx[j]) * delta_hh / xScaling[i];
386 }else{
387 DENSE_ELEM(Jac, j, i)((Jac->cols)[i][j]) = (fRes[j] - fx[j]) * delta_hh;
388 }
389
390 }
391 x[i] = xsave;
392 }
393
394 /* debug */
395 if (ACTIVE_STREAM(LOG_NLS_JAC)(useStream[LOG_NLS_JAC])){
396 infoStreamPrint(LOG_NLS_JAC, 1, "##KINSOL## omc dense matrix.");
397 PrintSparseMat(SlsConvertDls(Jac));
398 nlsKinsolJacSumDense(Jac);
399 messageClose(LOG_NLS_JAC);
400 }
401
402 /* performance measurement and statistics */
403 nlsData->jacobianTime += rt_ext_tp_tock(&(nlsData->jacobianTimeClock));
404 nlsData->numberOfJEval++;
405
406 return 0;
407}
408
409/* Element function for sparse matrix set */
410static void setJacElementKluSparse(int row, int col, double value, int nth, void* spJac)
411{
412 SlsMat mat = (SlsMat)spJac;
413 if (col > 0 && mat->colptrs[col] == 0){
414 mat->colptrs[col] = nth;
415 }
416 mat->rowvals[nth] = row;
417 mat->data[nth] = value;
418}
419
420/* finish sparse matrix, by fixing colprts */
421static void finishSparseColPtr(SlsMat mat)
422{
423 int i;
424 /* finish matrix colptrs */
425 mat->colptrs[mat->N] = mat->NNZ;
426 for(i=1; i<mat->N+1; ++i){
427 if (mat->colptrs[i] == 0){
428 warningStreamPrint(LOG_STDOUT, 0, "##KINSOL## Jacobian row %d singular. See LOG_NLS for more information.", i);
429 mat->colptrs[i] = mat->colptrs[i-1];
430 }
431 }
432}
433
434/*
435 * function calculates a jacobian matrix by
436 * numerical method finite differences with coloring
437 * into a sparse SlsMat matrix
438 */
439static
440int nlsSparseJac(N_Vector vecX, N_Vector vecFX, SlsMat Jac, void *userData, N_Vector tmp1, N_Vector tmp2)
441{
442 NLS_KINSOL_USERDATA *kinsolUserData = (NLS_KINSOL_USERDATA*) userData;
443 DATA* data = kinsolUserData->data;
444 threadData_t *threadData = kinsolUserData->threadData;
445 int sysNumber = kinsolUserData->sysNumber;
446 NONLINEAR_SYSTEM_DATA *nlsData = &(data->simulationInfo->nonlinearSystemData[sysNumber]);
447 NLS_KINSOL_DATA* kinsolData = (NLS_KINSOL_DATA*) nlsData->solverData;
448
449 /* prepare variables */
450 double *x = N_VGetArrayPointer(vecX);
451 double *fx = N_VGetArrayPointer(vecFX);
452 double *xsave = N_VGetArrayPointer(tmp1);
453 double *delta_hh = N_VGetArrayPointer(tmp2);
454 double *xScaling = NV_DATA_S(kinsolData->xScale)( ( (N_VectorContent_Serial)(kinsolData->xScale->content
) )->data )
;
455 double *fRes = NV_DATA_S(kinsolData->fRes)( ( (N_VectorContent_Serial)(kinsolData->fRes->content)
)->data )
;
456
457 SPARSE_PATTERN* sparsePattern = nlsData->sparsePattern;
458
459 const double delta_h = sqrt(DBL_EPSILON2.2204460492503131e-16*2e1);
460
461 long int i,j,ii;
462 int nth = 0;
463
464 /* performance measurement */
465 rt_ext_tp_tick(&nlsData->jacobianTimeClock);
466
467 /* reset matrix */
468 SlsSetToZero(Jac);
469
470 for(i = 0; i < sparsePattern->maxColors; i++)
471 {
472 for(ii=0; ii < kinsolData->size; ii++)
473 {
474 if(sparsePattern->colorCols[ii]-1 == i)
475 {
476 xsave[ii] = x[ii];
477 delta_hh[ii] = delta_h * (fabs(xsave[ii]) + 1.0);
478 if ((xsave[ii] + delta_hh[ii] >= nlsData->max[ii]))
479 delta_hh[ii] *= -1;
480 x[ii] += delta_hh[ii];
481
482 /* Calculate scaled difference quotient */
483 delta_hh[ii] = 1. / delta_hh[ii];
484 }
485 }
486 nlsKinsolResiduals(vecX, kinsolData->fRes, userData);
487
488 for(ii = 0; ii < kinsolData->size; ii++)
489 {
490 if(sparsePattern->colorCols[ii]-1 == i)
491 {
492 nth = sparsePattern->leadindex[ii];
493 while(nth < sparsePattern->leadindex[ii+1])
494 {
495 j = sparsePattern->index[nth];
496 if (kinsolData->nominalJac){
497 setJacElementKluSparse(j, ii, (fRes[j] - fx[j]) * delta_hh[ii] / xScaling[ii], nth, Jac);
498 }else{
499 setJacElementKluSparse(j, ii, (fRes[j] - fx[j]) * delta_hh[ii], nth, Jac);
500 }
501 nth++;
502 };
503 x[ii] = xsave[ii];
504 }
505 }
506 }
507 /* finish sparse matrix */
508 finishSparseColPtr(Jac);
509
510 /* debug */
511 if (ACTIVE_STREAM(LOG_NLS_JAC)(useStream[LOG_NLS_JAC])){
512 infoStreamPrint(LOG_NLS_JAC, 1, "##KINSOL## Sparse Matrix.");
513 PrintSparseMat(Jac);
514 nlsKinsolJacSumSparse(Jac);
515 messageClose(LOG_NLS_JAC);
516 }
517
518 /* performance measurement and statistics */
519 nlsData->jacobianTime += rt_ext_tp_tock(&(nlsData->jacobianTimeClock));
520 nlsData->numberOfJEval++;
521
522 return 0;
523}
524
525/*! \fn nlsSparseSymJac
526 *
527 * function calculates symbolic jacobian
528 *
529 * \param [ref] [data]
530 * \param [in] [sysNumber]
531
532 */
533static
534int nlsSparseSymJac(N_Vector vecX, N_Vector vecFX, SlsMat Jac, void *userData, N_Vector tmp1, N_Vector tmp2)
535{
536 NLS_KINSOL_USERDATA *kinsolUserData = (NLS_KINSOL_USERDATA*) userData;
537 DATA* data = kinsolUserData->data;
538 threadData_t *threadData = kinsolUserData->threadData;
539 int sysNumber = kinsolUserData->sysNumber;
540 NONLINEAR_SYSTEM_DATA *nlsData = &(data->simulationInfo->nonlinearSystemData[sysNumber]);
541 NLS_KINSOL_DATA* kinsolData = (NLS_KINSOL_DATA*) nlsData->solverData;
542
543 /* prepare variables */
544 double *x = N_VGetArrayPointer(vecX);
545 double *fx = N_VGetArrayPointer(vecFX);
546 double *xScaling = NV_DATA_S(kinsolData->xScale)( ( (N_VectorContent_Serial)(kinsolData->xScale->content
) )->data )
;
547
548 SPARSE_PATTERN* sparsePattern = nlsData->sparsePattern;
549 ANALYTIC_JACOBIAN* analyticJacobian = &data->simulationInfo->analyticJacobians[nlsData->jacobianIndex];
550
551 long int i,j,ii;
552 int nth = 0;
553 int nnz = sparsePattern->numberOfNoneZeros;
554
555 /* performance measurement */
556 rt_ext_tp_tick(&nlsData->jacobianTimeClock);
557
558 /* reset matrix */
559 SlsSetToZero(Jac);
560
561 if (analyticJacobian->constantEqns != NULL((void*)0)) {
562 analyticJacobian->constantEqns(data, threadData, analyticJacobian, NULL((void*)0));
563 }
564
565 for(i = 0; i < sparsePattern->maxColors; i++)
566 {
567 for(ii=0; ii < kinsolData->size; ii++)
568 {
569 if(sparsePattern->colorCols[ii]-1 == i)
570 {
571 analyticJacobian->seedVars[ii] = 1.0;
572 }
573 }
574 ((nlsData->analyticalJacobianColumn))(data, threadData, analyticJacobian, NULL((void*)0));
575
576 for(ii = 0; ii < kinsolData->size; ii++)
577 {
578 if(sparsePattern->colorCols[ii]-1 == i)
579 {
580 nth = sparsePattern->leadindex[ii];
581 while(nth < sparsePattern->leadindex[ii+1])
582 {
583 j = sparsePattern->index[nth];
584 if (kinsolData->nominalJac){
585 setJacElementKluSparse(j, ii, analyticJacobian->resultVars[j] / xScaling[ii] , nth, Jac);
586 }else{
587 setJacElementKluSparse(j, ii, analyticJacobian->resultVars[j], nth, Jac);
588 }
589 nth++;
590 };
591 analyticJacobian->seedVars[ii] = 0;
592 }
593 }
594 }
595 /* finish sparse matrix */
596 finishSparseColPtr(Jac);
597
598 /* debug */
599 if (ACTIVE_STREAM(LOG_NLS_JAC)(useStream[LOG_NLS_JAC])){
600 infoStreamPrint(LOG_NLS_JAC, 1, "##KINSOL## Sparse Matrix.");
601 PrintSparseMat(Jac);
602 nlsKinsolJacSumSparse(Jac);
603 messageClose(LOG_NLS_JAC);
604 }
605
606 /* performance measurement and statistics */
607 nlsData->jacobianTime += rt_ext_tp_tock(&(nlsData->jacobianTimeClock));
608 nlsData->numberOfJEval++;
609
610 return 0;
611}
612
613static
614void nlsKinsolJacSumDense(DlsMat mat)
615{
616 int i,j;
617 double sum;
618
619 for(i=0; i<mat->M; ++i){
620 sum = 0.0;
621 for(j=0; j<mat->N;++j){
622 sum += fabs(DENSE_ELEM(mat,j,i)((mat->cols)[i][j]));
623 }
624
625 if (sum == 0.0){
626 warningStreamPrint(LOG_NLS_V, 0, "sum of col %d of jacobian is zero!", i);
627 }else{
628 infoStreamPrint(LOG_NLS_JAC, 0, "col %d jac sum = %g", i, sum);
629 }
630
631 }
632}
633
634static
635void nlsKinsolJacSumSparse(SlsMat mat)
636{
637 int i,j;
638 double sum;
639
640 for(i=0; i<mat->N; ++i){
641 sum = 0.0;
642 for(j=mat->colptrs[i]; j<mat->colptrs[i+1];++j){
643 sum += fabs(mat->data[j]);
644 }
645
646 if (sum == 0.0){
647 warningStreamPrint(LOG_NLS_V, 0, "sum of col %d of jacobian is zero!", i);
648 }else{
649 infoStreamPrint(LOG_NLS_JAC, 0, "col %d jac sum = %g", i, sum);
650 }
651
652 }
653}
654
655
656static
657void nlsKinsolErrorPrint(int errorCode, const char *module, const char *function, char *msg, void *userData)
658{
659 NLS_KINSOL_DATA *kinsolData = (NLS_KINSOL_DATA*) userData;
660 DATA* data = kinsolData->userData.data;
661 int sysNumber = kinsolData->userData.sysNumber;
662 long eqSystemNumber = data->simulationInfo->nonlinearSystemData[sysNumber].equationIndex;
663
664 if (ACTIVE_STREAM(LOG_NLS_V)(useStream[LOG_NLS_V]))
665 {
666 warningStreamPrint(LOG_NLS_V, 1, "kinsol failed for %d", modelInfoGetEquation(&data->modelData->modelDataXml,eqSystemNumber).id);
667 warningStreamPrint(LOG_NLS_V, 0, "[module] %s | [function] %s | [error_code] %d", module, function, errorCode);
668 warningStreamPrint(LOG_NLS_V, 0, "%s", msg);
669
670 messageClose(LOG_NLS_V);
671 }
672}
673
674static
675void nlsKinsolInfoPrint(const char *module, const char *function, char *msg, void *user_data)
676{
677 if (ACTIVE_STREAM(LOG_NLS_V)(useStream[LOG_NLS_V]))
678 {
679 warningStreamPrint(LOG_NLS_V, 1, "%s %s:", module, function);
680 warningStreamPrint(LOG_NLS_V, 0, "%s", msg);
681
682 messageClose(LOG_NLS_V);
683 }
684}
685
686static
687void nlsKinsolSetMaxNewtonStep(NLS_KINSOL_DATA *kinsolData, double maxstepfactor)
688{
689 /* set maximum step size */
690 N_VConst(maxstepfactor, kinsolData->fTmp);
691 KINSetMaxNewtonStep(kinsolData->kinsolMemory, N_VWL2Norm(kinsolData->xScale, kinsolData->fTmp));
692
693}
694
695static
696void nlsKinsolResetInitial(DATA* data, NLS_KINSOL_DATA *kinsolData, NONLINEAR_SYSTEM_DATA* nlsData, int mode)
697{
698 double *xStart = NV_DATA_S(kinsolData->initialGuess)( ( (N_VectorContent_Serial)(kinsolData->initialGuess->
content) )->data )
;
699 /* set x vector */
700 switch (mode)
701 {
702 case INITIAL_EXTRAPOLATION0:
703 if(data->simulationInfo->discreteCall) {
704 memcpy(xStart, nlsData->nlsx, nlsData->size*(sizeof(double)));
705 } else {
706 memcpy(xStart, nlsData->nlsxExtrapolation, nlsData->size*(sizeof(double)));
707 }
708 break;
709 case INITIAL_OLDVALUES1:
710 memcpy(xStart, nlsData->nlsxOld, nlsData->size*(sizeof(double)));
711 break;
712 }
713}
714
715static
716void nlsKinsolXScaling(DATA* data, NLS_KINSOL_DATA *kinsolData, NONLINEAR_SYSTEM_DATA* nlsData, int mode)
717{
718 double *xStart = NV_DATA_S(kinsolData->initialGuess)( ( (N_VectorContent_Serial)(kinsolData->initialGuess->
content) )->data )
;
719 double *xScaling = NV_DATA_S(kinsolData->xScale)( ( (N_VectorContent_Serial)(kinsolData->xScale->content
) )->data )
;
720 int i;
721
722 /* if noScaling flag is used overwrite mode */
723 if (omc_flag[FLAG_NO_SCALING])
724 {
725 mode = SCALING_ONES2;
726 }
727
728 /* Use nominal value or the actual working point for scaling */
729 switch (mode)
730 {
731 case SCALING_NOMINALSTART1:
732 for (i=0;i<nlsData->size;i++){
733 xScaling[i] = 1.0/fmax(nlsData->nominal[i],fabs(xStart[i]));
734 }
735 break;
736 case SCALING_ONES2:
737 for (i=0;i<nlsData->size;i++){
738 xScaling[i] = 1.0;
739 }
740 break;
741 }
742}
743
744static
745void nlsKinsolFScaling(DATA* data, NLS_KINSOL_DATA *kinsolData, NONLINEAR_SYSTEM_DATA* nlsData, int mode)
746{
747 double *fScaling = NV_DATA_S(kinsolData->fScale)( ( (N_VectorContent_Serial)(kinsolData->fScale->content
) )->data )
;
748 N_Vector x = kinsolData->initialGuess;
749 int i,j;
750 SlsMat spJac;
751
752 /* if noScaling flag is used overwrite mode */
753 if (omc_flag[FLAG_NO_SCALING])
754 {
755 mode = SCALING_ONES2;
756 }
757
758 /* Use nominal value or the actual working point for scaling */
759 switch (mode)
760 {
761 case SCALING_JACOBIAN3:
762 {
763 N_Vector tmp1 = N_VNew_Serial(kinsolData->size);
764 N_Vector tmp2 = N_VNew_Serial(kinsolData->size);
765
766 /* enable scaled jacobian */
767 kinsolData->nominalJac = 1;
768
769 /* update x for the matrix */
770 nlsKinsolResiduals(x, kinsolData->fTmp, &kinsolData->userData);
771
772 /* calculate the right jacobian */
773 if(nlsData->isPatternAvailable && kinsolData->linearSolverMethod == NLS_LS_KLU)
774 {
775 spJac = NewSparseMat(kinsolData->size,kinsolData->size,kinsolData->nnz);
776 if (nlsData->analyticalJacobianColumn != NULL((void*)0)){
777 nlsSparseSymJac(x, kinsolData->fTmp, spJac, &kinsolData->userData, tmp1, tmp2);
778 } else {
779 nlsSparseJac(x, kinsolData->fTmp, spJac, &kinsolData->userData, tmp1, tmp2);
780 }
781 }
782 else
783 {
784 DlsMat denseJac = NewDenseMat(kinsolData->size,kinsolData->size);
785 nlsDenseJac(nlsData->size, x, kinsolData->fTmp, denseJac, &kinsolData->userData, tmp1, tmp2);
786 spJac = SlsConvertDls(denseJac);
787 }
788 /* disable scaled jacobian */
789 kinsolData->nominalJac = 0;
790
791 for (i=0;i<nlsData->size;i++){
792 fScaling[i] = 1e-12;
793 }
794 for(i=0; i<spJac->NNZ; ++i){
795 if (fScaling[spJac->rowvals[i]] < fabs(spJac->data[i])){
796 fScaling[spJac->rowvals[i]] = fabs(spJac->data[i]);
797 }
798 }
799 N_VInv(kinsolData->fScale, kinsolData->fScale);
800
801 DestroySparseMat(spJac);
802 N_VDestroy_Serial(tmp1);
803 N_VDestroy_Serial(tmp2);
804 break;
805 }
806 case SCALING_ONES2:
807 for (i=0;i<nlsData->size;i++){
808 fScaling[i] = 1.0;
809 }
810 break;
811 }
812}
813
814static
815void nlsKinsolConfigPrint(NLS_KINSOL_DATA *kinsolData, NONLINEAR_SYSTEM_DATA *nlsData)
816{
817 int retValue;
818 double fNorm;
819 DATA* data = kinsolData->userData.data;
820 int eqSystemNumber = nlsData->equationIndex;
821 _omc_vector vecStart, vecXScaling, vecFScaling;
822
823 if (!useStream[LOG_NLS_V]) {
824 return;
825 }
826
827 _omc_initVector(&vecStart, kinsolData->size, NV_DATA_S(kinsolData->initialGuess)( ( (N_VectorContent_Serial)(kinsolData->initialGuess->
content) )->data )
);
828 _omc_initVector(&vecXScaling, kinsolData->size, NV_DATA_S(kinsolData->xScale)( ( (N_VectorContent_Serial)(kinsolData->xScale->content
) )->data )
);
829 _omc_initVector(&vecFScaling, kinsolData->size, NV_DATA_S(kinsolData->fScale)( ( (N_VectorContent_Serial)(kinsolData->fScale->content
) )->data )
);
830
831 infoStreamPrint(LOG_NLS_V, 1, "Kinsol Configuration");
832 _omc_printVectorWithEquationInfo(&vecStart,
833 "Initial guess values", LOG_NLS_V, modelInfoGetEquation(&data->modelData->modelDataXml,eqSystemNumber));
834
835 _omc_printVectorWithEquationInfo(&vecXScaling,
836 "xScaling", LOG_NLS_V, modelInfoGetEquation(&data->modelData->modelDataXml,eqSystemNumber));
837
838 _omc_printVector(&vecFScaling, "fScaling", LOG_NLS_V);
839
840 infoStreamPrint(LOG_NLS_V, 0, "KINSOL F tolerance: %g", (*(KINMem)kinsolData->kinsolMemory).kin_fnormtol);
841 infoStreamPrint(LOG_NLS_V, 0, "KINSOL minimal step size %g", (*(KINMem)kinsolData->kinsolMemory).kin_scsteptol);
842 infoStreamPrint(LOG_NLS_V, 0, "KINSOL max iterations %d", 20*kinsolData->size);
843 infoStreamPrint(LOG_NLS_V, 0, "KINSOL strategy %d", kinsolData->kinsolStrategy);
844 infoStreamPrint(LOG_NLS_V, 0, "KINSOL current retry %d", kinsolData->retries);
845 infoStreamPrint(LOG_NLS_V, 0, "KINSOL max step %g", (*(KINMem)kinsolData->kinsolMemory).kin_mxnstepin);
846 infoStreamPrint(LOG_NLS_V, 0, "KINSOL linear solver %d", kinsolData->linearSolverMethod);
847
848 messageClose(LOG_NLS_V);
849}
850
851static
852int nlsKinsolErrorHandler(int errorCode, DATA *data, NONLINEAR_SYSTEM_DATA *nlsData, NLS_KINSOL_DATA *kinsolData)
853{
854 int retValue, i, retValue2=0, flag;
855 double fNorm;
856 double *xStart = NV_DATA_S(kinsolData->initialGuess)( ( (N_VectorContent_Serial)(kinsolData->initialGuess->
content) )->data )
;
857 double *xScaling = NV_DATA_S(kinsolData->xScale)( ( (N_VectorContent_Serial)(kinsolData->xScale->content
) )->data )
;
858 long outL;
859
860 /* check what kind of error
861 * retValue < 0 -> a non recoverable issue
862 * retValue == 1 -> try with other settings
863 */
864 KINSetNoInitSetup(kinsolData->kinsolMemory, FALSE0);
865
866 switch(errorCode)
867 {
868 case KIN_MEM_NULL-1:
869 case KIN_ILL_INPUT-2:
870 case KIN_NO_MALLOC-3:
871 errorStreamPrint(LOG_NLS_V, 0, "kinsol has a serious memory issue ERROR %d\n", errorCode);
872 return errorCode;
873 break;
874 /* just retry with new initial guess */
875 case KIN_MXNEWT_5X_EXCEEDED-7:
876 warningStreamPrint(LOG_NLS_V, 0, "Newton step exceed the maximum step size several times. Try again after increasing maximum step size.\n");
877 kinsolData->maxstepfactor *= 1e5;
878 nlsKinsolSetMaxNewtonStep(kinsolData, kinsolData->maxstepfactor);
879 return 1;
880 break;
881 /* just retry without line search */
882 case KIN_LINESEARCH_NONCONV-5:
883 warningStreamPrint(LOG_NLS_V, 0, "kinsols line search did not convergence. Try without.\n");
884 kinsolData->kinsolStrategy = KIN_NONE0;
885 kinsolData->retries--;
886 return 1;
887 /* maybe happened because of an out-dated factorization, so just retry */
888 case KIN_LSOLVE_FAIL-12:
889 warningStreamPrint(LOG_NLS_V, 0, "kinsols matrix need new factorization. Try again.\n");
890 if (nlsData->isPatternAvailable){
891 KINKLUReInit(kinsolData->kinsolMemory, kinsolData->size, kinsolData->nnz, 2);
892 }
893 return 1;
894 case KIN_MAXITER_REACHED-6:
895 case KIN_REPTD_SYSFUNC_ERR-15:
896 warningStreamPrint(LOG_NLS_V, 0, "kinsols runs into issues retry with different configuration.\n");
897 retValue = 1;
898 break;
899 case KIN_LSETUP_FAIL-11:
900 /* in case of something goes wrong with the symbolic jacobian try the numerical */
901 if ( kinsolData->linearSolverMethod == NLS_LS_KLU && nlsData->isPatternAvailable && nlsData->analyticalJacobianColumn != NULL((void*)0)){
902 flag = KINSlsSetSparseJacFn(kinsolData->kinsolMemory, nlsSparseJac);
903 }
904 if (checkReturnFlag(flag)){
905 errorStreamPrint(LOG_STDOUT, 0, "##KINSOL## Something goes wrong while initialize KINSOL Sparse Solver!");
906 return flag;
907 }
908 else
909 {
910 retValue = 1;
911 }
912 break;
913 case KIN_LINESEARCH_BCFAIL-8:
914 KINGetNumBetaCondFails(kinsolData->kinsolMemory, &outL);
915 warningStreamPrint(LOG_NLS_V, 0, "kinsols runs into issues with beta-condition fails: %ld\n", outL);
916 retValue = 1;
917 break;
918 default:
919 errorStreamPrint(LOG_STDOUT, 0, "kinsol has a serious solving issue ERROR %d\n", errorCode);
920 return errorCode;
921 break;
922 }
923
924 /* check if the current solution is sufficient anyway */
925 KINGetFuncNorm(kinsolData->kinsolMemory, &fNorm);
926 if (fNorm<FTOL_WITH_LESS_ACCURANCY1.e-6)
927 {
928 warningStreamPrint(LOG_NLS_V, 0, "Move forward with a less accurate solution.");
929 KINSetFuncNormTol(kinsolData->kinsolMemory, FTOL_WITH_LESS_ACCURANCY1.e-6);
930 KINSetScaledStepTol(kinsolData->kinsolMemory, FTOL_WITH_LESS_ACCURANCY1.e-6);
931 retValue2 = 1;
932 }
933 else
934 {
935 warningStreamPrint(LOG_NLS_V, 0, "Current status of fx = %f", fNorm);
936 }
937
938 /* reconfigure kinsol for an other try */
939 if (retValue == 1 && !retValue2)
940 {
941 switch(kinsolData->retries)
942 {
943 case 0:
944 /* try without scaling */
945 nlsKinsolXScaling(data, kinsolData, nlsData, SCALING_ONES2);
946 nlsKinsolFScaling(data, kinsolData, nlsData, SCALING_ONES2);
947 break;
948 case 1:
949 /* try without line-search and oldValues */
950 nlsKinsolResetInitial(data, kinsolData, nlsData, INITIAL_OLDVALUES1);
951 kinsolData->kinsolStrategy = KIN_LINESEARCH1;
952 break;
953 case 2:
954 /* try without line-search and oldValues */
955 nlsKinsolResetInitial(data, kinsolData, nlsData, INITIAL_EXTRAPOLATION0);
956 kinsolData->kinsolStrategy = KIN_NONE0;
957 break;
958 case 3:
959 /* try with exact newton */
960 nlsKinsolXScaling(data, kinsolData, nlsData, SCALING_NOMINALSTART1);
961 nlsKinsolFScaling(data, kinsolData, nlsData, SCALING_JACOBIAN3);
962 nlsKinsolResetInitial(data, kinsolData, nlsData, INITIAL_EXTRAPOLATION0);
963 KINSetMaxSetupCalls(kinsolData->kinsolMemory, 1);
964 kinsolData->kinsolStrategy = KIN_LINESEARCH1;
965 break;
966 case 4:
967 /* try with exact newton to with out x scaling values */
968 nlsKinsolXScaling(data, kinsolData, nlsData, SCALING_ONES2);
969 nlsKinsolFScaling(data, kinsolData, nlsData, SCALING_ONES2);
970 nlsKinsolResetInitial(data, kinsolData, nlsData, INITIAL_OLDVALUES1);
971 KINSetMaxSetupCalls(kinsolData->kinsolMemory, 1);
972 kinsolData->kinsolStrategy = KIN_LINESEARCH1;
973 break;
974 default:
975 retValue = 0;
976 break;
977 }
978 }
979
980 return retValue+retValue2;
981}
982
983int nlsKinsolSolve(DATA *data, threadData_t *threadData, int sysNumber)
984{
985 NONLINEAR_SYSTEM_DATA *nlsData = &(data->simulationInfo->nonlinearSystemData[sysNumber]);
986 NLS_KINSOL_DATA *kinsolData = (NLS_KINSOL_DATA*)nlsData->solverData;
987 long eqSystemNumber = nlsData->equationIndex;
988 int indexes[2] = {1,eqSystemNumber};
989
990 int flag, i;
991 long nFEval;
992 int success = 0;
993 int retry = 0;
994 double *xStart = NV_DATA_S(kinsolData->initialGuess)( ( (N_VectorContent_Serial)(kinsolData->initialGuess->
content) )->data )
;
995 double *xScaling = NV_DATA_S(kinsolData->xScale)( ( (N_VectorContent_Serial)(kinsolData->xScale->content
) )->data )
;
996 double *fScaling = NV_DATA_S(kinsolData->fScale)( ( (N_VectorContent_Serial)(kinsolData->fScale->content
) )->data )
;
997 double fNormValue;
998
999 /* set user data */
1000 kinsolData->userData.data = data;
1001 kinsolData->userData.threadData = threadData;
1002 kinsolData->userData.sysNumber = sysNumber;
1003
1004 infoStreamPrintWithEquationIndexes(LOG_NLS_V, 1, indexes, "Start Kinsol solver at time %g", data->localData[0]->timeValue);
1005
1006 /* */
1007 do{
1008 {
1009 /* It seems if we don't free KINSol on every iteration, it leaks memory.
1010 * But if we reset KINSol, it takes an enormous amount of time...
1011 */
1012 if (0) { /* Save time; leak memory */
1013 resetKinsolMemory(kinsolData, nlsData);
1014 } else {
1015 /* reset configuration settings */
1016 nlsKinsolConfigSetup(kinsolData);
1017 }
1018 }
1019
1020 nlsKinsolResetInitial(data, kinsolData, nlsData, INITIAL_EXTRAPOLATION0);
1021
1022 /* set x scaling */
1023 nlsKinsolXScaling(data, kinsolData, nlsData, SCALING_NOMINALSTART1);
1024
1025 /* set f scaling */
1026 nlsKinsolFScaling(data, kinsolData, nlsData, SCALING_JACOBIAN3);
1027
1028 /* set maximum step size */
1029 nlsKinsolSetMaxNewtonStep(kinsolData, kinsolData->maxstepfactor);
1030
1031 /* dump configuration */
1032 nlsKinsolConfigPrint(kinsolData, nlsData);
1033
1034 flag = KINSol(kinsolData->kinsolMemory, /* KINSol memory block */
1035 kinsolData->initialGuess, /* initial guess on input; solution vector */
1036 kinsolData->kinsolStrategy, /* global strategy choice */
1037 kinsolData->xScale, /* scaling vector, for the variable cc */
1038 kinsolData->fScale); /* scaling vector for function values fval */
1039
1040 infoStreamPrint(LOG_NLS_V, 0, "KINSol finished with errorCode %d.", flag);
1041 /* check for errors */
1042 if(flag < 0)
1043 {
1044 retry = nlsKinsolErrorHandler(flag, data, nlsData, kinsolData);
1045 }
1046
1047 /* solution found */
1048 if ((flag == KIN_SUCCESS0) || (flag == KIN_INITIAL_GUESS_OK1) || (flag == KIN_STEP_LT_STPTOL2))
1049 {
1050 success = 1;
1051 }
1052 kinsolData->retries++;
1053
1054 /* write statistics */
1055 KINGetNumNonlinSolvIters(kinsolData->kinsolMemory, &nFEval);
1056 nlsData->numberOfIterations += nFEval;
1057 nlsData->numberOfFEval += kinsolData->countResCalls;
1058
1059 infoStreamPrint(LOG_NLS_V, 0, "Next try? success = %d, retry = %d, retries = %d = %s\n", success, retry, kinsolData->retries,!success && !(retry<1) && kinsolData->retries<RETRY_MAX5 ? "true" : "false" );
1060 }while(!success && !(retry<0) && kinsolData->retries < RETRY_MAX5);
1061
1062 /* solution found */
1063 if (success)
1064 {
1065 /* check if solution really solves the residuals */
1066 nlsKinsolResiduals(kinsolData->initialGuess, kinsolData->fRes, &kinsolData->userData);
1067 if (!omc_flag[FLAG_NO_SCALING]){
1068 N_VProd(kinsolData->fRes, kinsolData->fScale, kinsolData->fRes);
1069 }
1070 fNormValue = N_VWL2Norm(kinsolData->fRes, kinsolData->fRes);
1071
1072 infoStreamPrint(LOG_NLS_V, 0, "%sEuclidean norm of F(u) = %e", (omc_flag[FLAG_NO_SCALING])?"":"scaled ", fNormValue);
1073 if (FTOL_WITH_LESS_ACCURANCY1.e-6<fNormValue)
1074 {
1075 warningStreamPrint(LOG_NLS_V, 0, "False positive solution. FNorm is not small enough.");
1076 success = 0;
1077 }
1078 else /* solved system for reuse linear solver information */
1079 {
1080 kinsolData->solved = 1;
1081 }
1082 /* copy solution */
1083 memcpy(nlsData->nlsx, xStart, nlsData->size*(sizeof(double)));
1084 }
1085
1086 messageClose(LOG_NLS_V);
1087
1088 return success;
1089}
1090
1091#else
1092
1093int nlsKinsolAllocate(int size, NONLINEAR_SYSTEM_DATA *nlsData, int linearSolverMethod)
1094{
1095 throwStreamPrint(NULL((void*)0),"no sundials/kinsol support activated");
1096 return 0;
1097}
1098
1099int nlsKinsolFree(void** solverData)
1100{
1101 throwStreamPrint(NULL((void*)0),"no sundials/kinsol support activated");
1102 return 0;
1103}
1104
1105int nlsKinsolSolve(DATA *data, threadData_t *threadData, int sysNumber)
1106{
1107 throwStreamPrint(threadData,"no sundials/kinsol support activated");
1108 return 0;
1109}
1110
1111#endif