1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
|
#include "alberta_util.h"
#include "alberta_util_intern.h"
/*--------------------------------------------------------------------------*/
/*--- information output for linear solvers ------------------------------*/
/*--------------------------------------------------------------------------*/
void start_info(const char *funcName, OEM_DATA *oem)
{
oem->info = oem->info > 10 ? 10 : oem->info;
INFO(oem->info,6,"with tolerance %le", oem->tolerance);
if (oem->restart > 0)
PRINT_INFO(oem->info,6," and restart %d\n", oem->restart);
else
PRINT_INFO(oem->info,6,"\n");
INFO(oem->info,2,"iter. | residual | red.\n");
fflush(stdout);
return;
}
void break_info(const char *funcName, OEM_DATA *oem, const char *reason,
int iter, REAL res, REAL *ores, WORKSPACE *ws)
{
if (*ores && *ores > 0)
INFO(oem->info,2,"%5d | %12.5le | %8.2le\n", iter, res, res/(*ores));
else
INFO(oem->info,2,"%5d | %12.5le |\n", iter);
INFO(oem->info,2,"stop due to: %s\n", reason);
fflush(stdout);
free_oem_workspace(ws, oem);
oem->residual = res;
}
int solve_info(const char *funcName, OEM_DATA *oem, int iter, REAL res,
REAL *ores, WORKSPACE *ws)
{
static int step[11] = {0, 1000, 500, 200, 100, 50, 20, 10, 5, 2, 1};
if (res <= oem->tolerance || (oem->info && (iter%step[oem->info] == 0))
|| iter == oem->max_iter)
{
if (*ores)
{
if (*ores > 0.0)
{
REAL red = res/(*ores);
INFO(oem->info,2,"%5d | %12.5le | %8.2le\n", iter, res, red);
}
else
{
INFO(oem->info,2,"%5d | %12.5le | --------\n", iter, res);
}
*ores = res;
}
else
{
INFO(oem->info,2,"%5d | %12.5le |\n", iter, res);
}
}
oem->residual = res;
if (iter >= oem->max_iter || res <= oem->tolerance) {
if (res > oem->tolerance) {
INFO(oem->info,1,"tolerance %le not reached after %d iterations\n",
oem->tolerance, iter);
} else {
INFO(oem->info,6,"finished successfully with %d iterations\n",iter);
}
fflush(stdout);
free_oem_workspace(ws, oem);
return 1;
}
fflush(stdout);
return 0;
}
/*--------------------------------------------------------------------------*/
/*--- checking of workspace, reallocation of workspace if neccesary ------*/
/*--------------------------------------------------------------------------*/
WORKSPACE *check_workspace(const char *funcName, const char *file, int line,
size_t size, WORKSPACE *ws)
{
if (!ws)
{
ws = GET_WORKSPACE(size*sizeof(REAL));
}
else if (size*sizeof(REAL) > ws->size)
{
WARNING("need workspace for %d REALs\n", size);
WARNING("reallocating workspace of length %d\n", size*sizeof(REAL));
REALLOC_WORKSPACE(ws, size*sizeof(REAL));
}
return(ws);
}
void free_oem_workspace(WORKSPACE *ws, OEM_DATA *oem)
{
if (ws != oem->ws)
FREE_WORKSPACE(ws);
return;
}
|
