source: trunk/Clp/src/ClpHelperFunctions.cpp @ 1722

Last change on this file since 1722 was 1722, checked in by stefan, 8 years ago

adjust to changes in CoinUtils? header files

  • Property svn:eol-style set to native
  • Property svn:keywords set to Id
File size: 8.0 KB
Line 
1/* $Id: ClpHelperFunctions.cpp 1722 2011-04-17 09:58:37Z stefan $ */
2// Copyright (C) 2003, International Business Machines
3// Corporation and others.  All Rights Reserved.
4// This code is licensed under the terms of the Eclipse Public License (EPL).
5
6/*
7    Note (JJF) I have added some operations on arrays even though they may
8    duplicate CoinDenseVector.  I think the use of templates was a mistake
9    as I don't think inline generic code can take as much advantage of
10    parallelism or machine architectures or memory hierarchies.
11
12*/
13#include <cfloat>
14#include <cstdlib>
15#include <cmath>
16#include "CoinHelperFunctions.hpp"
17#include "CoinTypes.hpp"
18
19double
20maximumAbsElement(const double * region, int size)
21{
22     int i;
23     double maxValue = 0.0;
24     for (i = 0; i < size; i++)
25          maxValue = CoinMax(maxValue, fabs(region[i]));
26     return maxValue;
27}
28void
29setElements(double * region, int size, double value)
30{
31     int i;
32     for (i = 0; i < size; i++)
33          region[i] = value;
34}
35void
36multiplyAdd(const double * region1, int size, double multiplier1,
37            double * region2, double multiplier2)
38{
39     int i;
40     if (multiplier1 == 1.0) {
41          if (multiplier2 == 1.0) {
42               for (i = 0; i < size; i++)
43                    region2[i] = region1[i] + region2[i];
44          } else if (multiplier2 == -1.0) {
45               for (i = 0; i < size; i++)
46                    region2[i] = region1[i] - region2[i];
47          } else if (multiplier2 == 0.0) {
48               for (i = 0; i < size; i++)
49                    region2[i] = region1[i] ;
50          } else {
51               for (i = 0; i < size; i++)
52                    region2[i] = region1[i] + multiplier2 * region2[i];
53          }
54     } else if (multiplier1 == -1.0) {
55          if (multiplier2 == 1.0) {
56               for (i = 0; i < size; i++)
57                    region2[i] = -region1[i] + region2[i];
58          } else if (multiplier2 == -1.0) {
59               for (i = 0; i < size; i++)
60                    region2[i] = -region1[i] - region2[i];
61          } else if (multiplier2 == 0.0) {
62               for (i = 0; i < size; i++)
63                    region2[i] = -region1[i] ;
64          } else {
65               for (i = 0; i < size; i++)
66                    region2[i] = -region1[i] + multiplier2 * region2[i];
67          }
68     } else if (multiplier1 == 0.0) {
69          if (multiplier2 == 1.0) {
70               // nothing to do
71          } else if (multiplier2 == -1.0) {
72               for (i = 0; i < size; i++)
73                    region2[i] =  -region2[i];
74          } else if (multiplier2 == 0.0) {
75               for (i = 0; i < size; i++)
76                    region2[i] =  0.0;
77          } else {
78               for (i = 0; i < size; i++)
79                    region2[i] =  multiplier2 * region2[i];
80          }
81     } else {
82          if (multiplier2 == 1.0) {
83               for (i = 0; i < size; i++)
84                    region2[i] = multiplier1 * region1[i] + region2[i];
85          } else if (multiplier2 == -1.0) {
86               for (i = 0; i < size; i++)
87                    region2[i] = multiplier1 * region1[i] - region2[i];
88          } else if (multiplier2 == 0.0) {
89               for (i = 0; i < size; i++)
90                    region2[i] = multiplier1 * region1[i] ;
91          } else {
92               for (i = 0; i < size; i++)
93                    region2[i] = multiplier1 * region1[i] + multiplier2 * region2[i];
94          }
95     }
96}
97double
98innerProduct(const double * region1, int size, const double * region2)
99{
100     int i;
101     double value = 0.0;
102     for (i = 0; i < size; i++)
103          value += region1[i] * region2[i];
104     return value;
105}
106void
107getNorms(const double * region, int size, double & norm1, double & norm2)
108{
109     norm1 = 0.0;
110     norm2 = 0.0;
111     int i;
112     for (i = 0; i < size; i++) {
113          norm2 += region[i] * region[i];
114          norm1 = CoinMax(norm1, fabs(region[i]));
115     }
116}
117#if COIN_LONG_WORK
118// For long double versions
119CoinWorkDouble
120maximumAbsElement(const CoinWorkDouble * region, int size)
121{
122     int i;
123     CoinWorkDouble maxValue = 0.0;
124     for (i = 0; i < size; i++)
125          maxValue = CoinMax(maxValue, CoinAbs(region[i]));
126     return maxValue;
127}
128void
129setElements(CoinWorkDouble * region, int size, CoinWorkDouble value)
130{
131     int i;
132     for (i = 0; i < size; i++)
133          region[i] = value;
134}
135void
136multiplyAdd(const CoinWorkDouble * region1, int size, CoinWorkDouble multiplier1,
137            CoinWorkDouble * region2, CoinWorkDouble multiplier2)
138{
139     int i;
140     if (multiplier1 == 1.0) {
141          if (multiplier2 == 1.0) {
142               for (i = 0; i < size; i++)
143                    region2[i] = region1[i] + region2[i];
144          } else if (multiplier2 == -1.0) {
145               for (i = 0; i < size; i++)
146                    region2[i] = region1[i] - region2[i];
147          } else if (multiplier2 == 0.0) {
148               for (i = 0; i < size; i++)
149                    region2[i] = region1[i] ;
150          } else {
151               for (i = 0; i < size; i++)
152                    region2[i] = region1[i] + multiplier2 * region2[i];
153          }
154     } else if (multiplier1 == -1.0) {
155          if (multiplier2 == 1.0) {
156               for (i = 0; i < size; i++)
157                    region2[i] = -region1[i] + region2[i];
158          } else if (multiplier2 == -1.0) {
159               for (i = 0; i < size; i++)
160                    region2[i] = -region1[i] - region2[i];
161          } else if (multiplier2 == 0.0) {
162               for (i = 0; i < size; i++)
163                    region2[i] = -region1[i] ;
164          } else {
165               for (i = 0; i < size; i++)
166                    region2[i] = -region1[i] + multiplier2 * region2[i];
167          }
168     } else if (multiplier1 == 0.0) {
169          if (multiplier2 == 1.0) {
170               // nothing to do
171          } else if (multiplier2 == -1.0) {
172               for (i = 0; i < size; i++)
173                    region2[i] =  -region2[i];
174          } else if (multiplier2 == 0.0) {
175               for (i = 0; i < size; i++)
176                    region2[i] =  0.0;
177          } else {
178               for (i = 0; i < size; i++)
179                    region2[i] =  multiplier2 * region2[i];
180          }
181     } else {
182          if (multiplier2 == 1.0) {
183               for (i = 0; i < size; i++)
184                    region2[i] = multiplier1 * region1[i] + region2[i];
185          } else if (multiplier2 == -1.0) {
186               for (i = 0; i < size; i++)
187                    region2[i] = multiplier1 * region1[i] - region2[i];
188          } else if (multiplier2 == 0.0) {
189               for (i = 0; i < size; i++)
190                    region2[i] = multiplier1 * region1[i] ;
191          } else {
192               for (i = 0; i < size; i++)
193                    region2[i] = multiplier1 * region1[i] + multiplier2 * region2[i];
194          }
195     }
196}
197CoinWorkDouble
198innerProduct(const CoinWorkDouble * region1, int size, const CoinWorkDouble * region2)
199{
200     int i;
201     CoinWorkDouble value = 0.0;
202     for (i = 0; i < size; i++)
203          value += region1[i] * region2[i];
204     return value;
205}
206void
207getNorms(const CoinWorkDouble * region, int size, CoinWorkDouble & norm1, CoinWorkDouble & norm2)
208{
209     norm1 = 0.0;
210     norm2 = 0.0;
211     int i;
212     for (i = 0; i < size; i++) {
213          norm2 += region[i] * region[i];
214          norm1 = CoinMax(norm1, CoinAbs(region[i]));
215     }
216}
217#endif
218#ifdef DEBUG_MEMORY
219#include <malloc.h>
220#include <stdio.h>
221#include <stdlib.h>
222
223typedef void (*NEW_HANDLER)();
224static NEW_HANDLER new_handler;                        // function to call if `new' fails (cf. ARM p. 281)
225
226// Allocate storage.
227void *
228operator new(size_t size)
229{
230     void * p;
231     for (;;) {
232          p = malloc(size);
233          if      (p)           break;        // success
234          else if (new_handler) new_handler();   // failure - try again (allow user to release some storage first)
235          else                  break;        // failure - no retry
236     }
237     if (size > 1000000)
238          printf("Allocating memory of size %d\n", size);
239     return p;
240}
241
242// Deallocate storage.
243void
244operator delete(void *p)
245{
246     free(p);
247     return;
248}
249void
250operator delete [] (void *p)
251{
252     free(p);
253     return;
254}
255#endif
Note: See TracBrowser for help on using the repository browser.