Actual source code: tvecsingle.c

  1: #include "tvecsingle.h"
  2: #include "tao_general.h"
  3: #include "stdio.h"

  5: TaoVecFloatArray::TaoVecFloatArray(TaoInt nn):TaoVec(),n(nn){
  6:   v=new float[nn];
  7:   fallocated=1;
  8:   return;
  9: }

 11: TaoVecFloatArray::TaoVecFloatArray(TaoInt nn, float *ff):TaoVec(),n(nn){
 12:   v=ff;
 13:   fallocated=0;
 14:   return;
 15: }

 17: int TaoVecFloatArray::Clone( TaoVec** tv ){

 19:   *tv = new TaoVecFloatArray(this->n);
 20:   int info = (*tv)->CopyFrom(this);CHKERRQ(info);
 21:   return 0;
 22: }

 24: int TaoVecFloatArray::GetFloats(float **fptr, TaoInt *nn){
 25:   *fptr=v;
 26:   *nn=n;
 27:   return 0;
 28: }

 30: int TaoVecFloatArray::RestoreFloats(float **fptr, TaoInt *nn){
 31:   *fptr=0;
 32:   *nn=0;
 33:   return 0;
 34: }

 36: int TaoVecFloatArray::GetDimension(TaoInt *nn){
 37:   *nn=n;
 38:   return 0;
 39: }

 41: int TaoVecFloatArray::Compatible(TaoVec *tv, TaoTruth *flag){
 42:   TaoInt nn;
 43:   int info;
 44:   float *fptr;
 45:   TaoVecFloatArray* vv =  (TaoVecFloatArray*)(tv);

 47:   info = vv->GetData(&fptr,&nn);
 48:   if (info==0 && nn == n) *flag=TAO_TRUE;
 49:   else *flag=TAO_FALSE;
 50:   return 0;
 51: }

 53: int TaoVecFloatArray::View(){
 54:   for (TaoInt i=0;i<n;++i)
 55:     printf(" %4.2e \n ",v[i]);
 56:   return 0;
 57: }


 62: int TaoVecFloatArray::SetToConstant( double c ){
 63:   TaoInt i,nn;
 64:   int info;
 65:   float *tptr;

 67:   TaoFunctionBegin;
 68:   info = this->GetFloats(&tptr,&nn);CHKERRQ(info);
 69:   for (i=0;i<nn;i++){ tptr[i]=c; }
 70:   info = this->RestoreFloats(&tptr,&nn);CHKERRQ(info);
 71:   TaoFunctionReturn(0);
 72: }

 76: int TaoVecFloatArray::SetToZero(){
 77:   TaoInt i,nn;
 78:   int info;
 79:   float *tptr;

 81:   TaoFunctionBegin;
 82:   info = this->GetFloats(&tptr,&nn);CHKERRQ(info);
 83:   for (i=0;i<nn;i++){ tptr[i]=0; }
 84:   info = this->RestoreFloats(&tptr,&nn);CHKERRQ(info);
 85:   TaoFunctionReturn(0);
 86: }

 90: int TaoVecFloatArray::CopyFrom( TaoVec* tv ){
 91:   TaoInt i,nn1,nn2;
 92:   int info;
 93:   float *tptr1,*tptr2;
 94:   TaoVecFloatArray* vv =  (TaoVecFloatArray*)(tv);

 96:   TaoFunctionBegin;
 97:   info = this->GetFloats(&tptr1,&nn1);CHKERRQ(info);
 98:   info = vv->GetFloats(&tptr2,&nn2);CHKERRQ(info);

100:   for (i=0;i<nn1;i++){ tptr1[i]=tptr2[i]; }

102:   info = vv->RestoreFloats(&tptr2,&nn2);CHKERRQ(info);
103:   info = this->RestoreFloats(&tptr1,&nn1);CHKERRQ(info);
104:   TaoFunctionReturn(0);
105: }

109: int TaoVecFloatArray::ScaleCopyFrom( double a, TaoVec* tv ){
110:   TaoInt i,nn1,nn2;
111:   int info;
112:   float *tptr1,*tptr2;
113:   TaoVecFloatArray* vv =  (TaoVecFloatArray*)(tv);

115:   TaoFunctionBegin;
116:   info = this->GetFloats(&tptr1,&nn1);CHKERRQ(info);
117:   info = vv->GetFloats(&tptr2,&nn2);CHKERRQ(info);
118:   for (i=0;i<nn1;i++){ tptr1[i]=a*tptr2[i]; }
119:   info = vv->RestoreFloats(&tptr2,&nn2);CHKERRQ(info);
120:   info = this->RestoreFloats(&tptr1,&nn1);CHKERRQ(info);
121:   TaoFunctionReturn(0);
122: }

126: int TaoVecFloatArray::NormInfinity(double *vnorm){
127:   TaoInt i,nn;
128:   int info;
129:   float dd=0, *vv;

131:   TaoFunctionBegin;
132:   info = this->GetFloats(&vv,&nn);CHKERRQ(info);
133:   for (i=0;i<nn;i++){
134:     if (vv[i]<0 && vv[i]<-dd) dd=-vv[i];
135:     else if (vv[i]>0 && vv[i]>dd) dd=vv[i];
136:   }
137:   info = this->RestoreFloats(&vv,&nn);CHKERRQ(info);
138:   *vnorm=dd;
139:   TaoFunctionReturn(0);
140: }

144: int TaoVecFloatArray::Norm1(double *vnorm){
145:   TaoInt i,nn;
146:   int info;
147:   float dd=0, *vv;

149:   TaoFunctionBegin;
150:   info = this->GetFloats(&vv,&nn);CHKERRQ(info);
151:   for (i=0;i<nn;i++){ 
152:     if (vv[i]<0) dd-=vv[i]; 
153:     else dd+=vv[i];
154:   }
155:   info = this->RestoreFloats(&vv,&nn);CHKERRQ(info);
156:   TaoFunctionReturn(0);
157: }

161: int TaoVecFloatArray::Norm2(double *vnorm){
162:   TaoInt i,nn;
163:   int info;
164:   float dd=0, *vv;

166:   TaoFunctionBegin;
167:   info = this->GetFloats(&vv,&nn);CHKERRQ(info);
168:   for (i=0;i<nn;i++) dd+=vv[i]*vv[i];
169:   *vnorm=sqrt(dd);
170:   info = this->RestoreFloats(&vv,&nn);CHKERRQ(info);
171:   TaoFunctionReturn(0);
172: }

176: int TaoVecFloatArray::Norm2squared(double *vnorm2){
177:   TaoInt i,nn;
178:   int info;
179:   float dd=0, *vv;

181:   TaoFunctionBegin;
182:   info = this->GetFloats(&vv,&nn);CHKERRQ(info);
183:   for (i=0;i<nn;i++) dd+=vv[i]*vv[i];
184:   *vnorm2=dd;
185:   info = this->RestoreFloats(&vv,&nn);CHKERRQ(info);
186:   TaoFunctionReturn(0);
187: }

191: int TaoVecFloatArray::Scale( double alpha ){
192:   TaoInt i,nn;
193:   int info;
194:   float *vv;

196:   TaoFunctionBegin;
197:   info = this->GetFloats(&vv,&nn);CHKERRQ(info);
198:   for (i=0;i<nn;i++) vv[i]*=alpha;
199:   info = this->RestoreFloats(&vv,&nn);CHKERRQ(info);
200:   TaoFunctionReturn(0);
201: }

205: int TaoVecFloatArray::Axpy( double alpha, TaoVec* tv ){
206:   TaoInt i,nn1,nn2;
207:   int info;
208:   float *tptr1,*tptr2;
209:   TaoVecFloatArray* xx =  (TaoVecFloatArray*)(tv);

211:   TaoFunctionBegin;
212:   info = this->GetFloats(&tptr1,&nn1);CHKERRQ(info);
213:   info = xx->GetFloats(&tptr2,&nn2);CHKERRQ(info);
214:   for (i=0;i<nn1;i++){ tptr1[i]+= alpha * tptr2[i]; }
215:   info = this->RestoreFloats(&tptr1,&nn1);CHKERRQ(info);
216:   info = xx->RestoreFloats(&tptr2,&nn2);CHKERRQ(info);
217:   TaoFunctionReturn(0);
218: }

222: int TaoVecFloatArray::Axpby( double alpha, TaoVec* tv, double beta ){
223:   TaoInt i,nn1,nn2;
224:   int info;
225:   float *tptr1,*tptr2;
226:   TaoVecFloatArray* xx =  (TaoVecFloatArray*)(tv);

228:   TaoFunctionBegin;
229:   info = this->GetFloats(&tptr1,&nn1);CHKERRQ(info);
230:   info = xx->GetFloats(&tptr2,&nn2);CHKERRQ(info);
231:   for (i=0;i<nn1;i++){ tptr1[i] = beta * tptr1[i] + alpha * tptr2[i]; }
232:   info = this->RestoreFloats(&tptr1,&nn1);CHKERRQ(info);
233:   info = xx->RestoreFloats(&tptr2,&nn2);CHKERRQ(info);
234:   TaoFunctionReturn(0);
235: }

239: int TaoVecFloatArray::Aypx( double alpha, TaoVec* tv ){
240:   TaoInt i,nn1,nn2;
241:   int info;
242:   float *tptr1,*tptr2;
243:   TaoVecFloatArray* xx =  (TaoVecFloatArray*)(tv);

245:   TaoFunctionBegin;
246:   info = this->GetFloats(&tptr1,&nn1);CHKERRQ(info);
247:   info = xx->GetFloats(&tptr2,&nn2);CHKERRQ(info);
248:   for (i=0;i<nn1;i++){ tptr1[i] = alpha * tptr1[i] + tptr2[i]; }
249:   info = this->RestoreFloats(&tptr1,&nn1);CHKERRQ(info);
250:   info = xx->RestoreFloats(&tptr2,&nn2);CHKERRQ(info);
251:   TaoFunctionReturn(0);
252: }
253:  
256: int TaoVecFloatArray::AddConstant( double alpha ){
257:   TaoInt i,nn;
258:   int info;
259:   float *vv;

261:   TaoFunctionBegin;
262:   info = this->GetFloats(&vv,&nn);CHKERRQ(info);
263:   for (i=0;i<nn;i++) vv[i]+=alpha;
264:   info = this->RestoreFloats(&vv,&nn);CHKERRQ(info);
265:   TaoFunctionReturn(0);
266: }

270: int TaoVecFloatArray::Dot( TaoVec* tv, double *vDotv ){
271:   TaoInt i,nn1,nn2;
272:   int info;
273:   float dd=0,*tptr1,*tptr2;
274:   TaoVecFloatArray* vv =  (TaoVecFloatArray*)(tv);

276:   TaoFunctionBegin;
277:   info = this->GetFloats(&tptr1,&nn1);CHKERRQ(info);
278:   info = vv->GetFloats(&tptr2,&nn2);CHKERRQ(info);
279:   for (i=0;i<nn1;i++) dd+=tptr1[i]*tptr2[i];
280:   info = this->RestoreFloats(&tptr1,&nn1);CHKERRQ(info);
281:   info = vv->RestoreFloats(&tptr2,&nn2);CHKERRQ(info);
282:   *vDotv=dd;
283:   TaoFunctionReturn(0);
284: }

288: int TaoVecFloatArray::Negate(){ 
289:   TaoInt i,nn;
290:   int info;
291:   float *vv;

293:   TaoFunctionBegin;
294:   info = this->GetFloats(&vv,&nn);CHKERRQ(info);
295:   for (i=0;i<nn;i++) vv[i]=-vv[i];
296:   info = this->RestoreFloats(&vv,&nn);CHKERRQ(info);
297:   TaoFunctionReturn(0);
298: }

302: int TaoVecFloatArray::Reciprocal(){ 
303:   TaoInt i,nn;
304:   int info;
305:   float *vv;

307:   TaoFunctionBegin;
308:   info = this->GetFloats(&vv,&nn);CHKERRQ(info);
309:   for (i=0;i<nn;i++){ 
310:     if (vv[i]!=0)  vv[i]= 1.0/vv[i];
311:     else vv[i]=TAO_INFINITY;
312:   }
313:   info = this->RestoreFloats(&vv,&nn);CHKERRQ(info);
314:   TaoFunctionReturn(0);
315: }

319: int TaoVecFloatArray::Sqrt(){ 
320:   TaoInt i,nn;
321:   int info;
322:   float *vv;

324:   TaoFunctionBegin;
325:   info = this->GetFloats(&vv,&nn);CHKERRQ(info);
326:   for (i=0;i<nn;i++){ 
327:     if (vv[i] >= 0)  vv[i]= sqrt(vv[i]);
328:     else vv[i]=TAO_INFINITY;
329:   }
330:   info = this->RestoreFloats(&vv,&nn);CHKERRQ(info);
331:   TaoFunctionReturn(0);
332: }

336: int TaoVecFloatArray::Pow(double p){ 
337:   TaoInt i,nn;
338:   int info;
339:   float *vv;

341:   TaoFunctionBegin;
342:   info = this->GetFloats(&vv,&nn);CHKERRQ(info);
343:   for (i=0;i<nn;i++){ 
344:     if (vv[i] >= 0)  vv[i]= pow((float)vv[i], (float)p);
345:     else vv[i]=TAO_INFINITY;
346:   }
347:   info = this->RestoreFloats(&vv,&nn);CHKERRQ(info);
348:   TaoFunctionReturn(0);
349: }

353: int TaoVecFloatArray::PointwiseMultiply( TaoVec* tv, TaoVec* tw ){
354:   TaoInt i,nn1,nn2,nn3;
355:   int info;
356:   float *tptr1,*tptr2,*tptr3;
357:   TaoVecFloatArray* vv =  (TaoVecFloatArray*)(tv);
358:   TaoVecFloatArray* ww =  (TaoVecFloatArray*)(tw);

360:   TaoFunctionBegin;
361:   info = this->GetFloats(&tptr1,&nn1);CHKERRQ(info);
362:   info = vv->GetFloats(&tptr2,&nn2);CHKERRQ(info);
363:   info = ww->GetFloats(&tptr3,&nn3);CHKERRQ(info);
364:   for (i=0;i<nn1;i++) tptr1[i]=tptr2[i] * tptr3[i];
365:   info = this->RestoreFloats(&tptr1,&nn1);CHKERRQ(info);
366:   info = vv->RestoreFloats(&tptr2,&nn2);CHKERRQ(info);
367:   info = ww->RestoreFloats(&tptr3,&nn3);CHKERRQ(info);

369:   TaoFunctionReturn(0);
370: }

374: int TaoVecFloatArray::PointwiseDivide( TaoVec* tv , TaoVec* tw){
375:   TaoInt i,nn1,nn2,nn3;
376:   int info;
377:   float *tptr1,*tptr2,*tptr3;
378:   TaoVecFloatArray* vv =  (TaoVecFloatArray*)(tv);
379:   TaoVecFloatArray* ww =  (TaoVecFloatArray*)(tw);

381:   TaoFunctionBegin;
382:   info = this->GetFloats(&tptr1,&nn1);CHKERRQ(info);
383:   info = vv->GetFloats(&tptr2,&nn2);CHKERRQ(info);
384:   info = ww->GetFloats(&tptr3,&nn3);CHKERRQ(info);

386:   for (i=0;i<nn1;i++) tptr1[i]=tptr2[i] / tptr3[i];
387:   info = this->RestoreFloats(&tptr1,&nn1);CHKERRQ(info);
388:   info = vv->RestoreFloats(&tptr2,&nn2);CHKERRQ(info);
389:   info = ww->RestoreFloats(&tptr3,&nn3);CHKERRQ(info);

391:   TaoFunctionReturn(0);
392: }

396: int TaoVecFloatArray::Median( TaoVec* tv, TaoVec* tw, TaoVec* tx){
397:   TaoInt i,nn1,nn2,nn3,nn4;
398:   int info;
399:   float *tptr1,*tptr2,*tptr3,*tptr4;
400:   TaoVecFloatArray* vv =  (TaoVecFloatArray*)(tv);
401:   TaoVecFloatArray* ww =  (TaoVecFloatArray*)(tw);
402:   TaoVecFloatArray* xx =  (TaoVecFloatArray*)(tx);

404:   TaoFunctionBegin;
405:   info = this->GetFloats(&tptr1,&nn1);CHKERRQ(info);
406:   info = vv->GetFloats(&tptr2,&nn2);CHKERRQ(info);
407:   info = ww->GetFloats(&tptr3,&nn3);CHKERRQ(info);
408:   info = xx->GetFloats(&tptr4,&nn4);CHKERRQ(info);

410:   for (i=0;i<nn1;i++){
411:     if (tptr2[i]<=tptr3[i] && tptr3[i] <= tptr4[i]){
412:       tptr1[i]=tptr3[i];
413:     } else if (tptr4[i]<=tptr3[i] && tptr3[i] <= tptr2[i]){
414:       tptr1[i]=tptr3[i];
415:     } else if (tptr3[i]<=tptr2[i] && tptr2[i] <= tptr4[i]){
416:       tptr1[i]=tptr2[i];
417:     } else if (tptr4[i]<=tptr2[i] && tptr2[i] <= tptr3[i]){
418:       tptr1[i]=tptr2[i];
419:     } else {
420:       tptr1[i]=tptr4[i];
421:     }
422:   }
423:   info = this->RestoreFloats(&tptr1,&nn1);CHKERRQ(info);
424:   info = vv->RestoreFloats(&tptr2,&nn2);CHKERRQ(info);
425:   info = ww->RestoreFloats(&tptr3,&nn3);CHKERRQ(info);
426:   info = xx->RestoreFloats(&tptr4,&nn4);CHKERRQ(info);

428:   TaoFunctionReturn(0);
429: }

433: int TaoVecFloatArray::PointwiseMinimum( TaoVec* tv, TaoVec* tw){
434:   TaoInt i,nn1,nn2,nn3;
435:   int info;
436:   float *tptr1,*tptr2,*tptr3;
437:   TaoVecFloatArray* vv =  (TaoVecFloatArray*)(tv);
438:   TaoVecFloatArray* ww =  (TaoVecFloatArray*)(tw);

440:   TaoFunctionBegin;
441:   info = this->GetFloats(&tptr1,&nn1);CHKERRQ(info);
442:   info = vv->GetFloats(&tptr2,&nn2);CHKERRQ(info);
443:   info = ww->GetFloats(&tptr3,&nn3);CHKERRQ(info);

445:   for (i=0;i<nn1;i++) tptr1[i] = TaoMin( tptr2[i] , tptr3[i]);

447:   info = this->RestoreFloats(&tptr1,&nn1);CHKERRQ(info);
448:   info = vv->RestoreFloats(&tptr2,&nn2);CHKERRQ(info);
449:   info = ww->RestoreFloats(&tptr3,&nn3);CHKERRQ(info);

451:   TaoFunctionReturn(0);
452: }

456: int TaoVecFloatArray::PointwiseMaximum( TaoVec* tv, TaoVec* tw){
457:   TaoInt i,nn1,nn2,nn3;
458:   int info;
459:   float *tptr1,*tptr2,*tptr3;
460:   TaoVecFloatArray* vv =  (TaoVecFloatArray*)(tv);
461:   TaoVecFloatArray* ww =  (TaoVecFloatArray*)(tw);

463:   TaoFunctionBegin;
464:   info = this->GetFloats(&tptr1,&nn1);CHKERRQ(info);
465:   info = vv->GetFloats(&tptr2,&nn2);CHKERRQ(info);
466:   info = ww->GetFloats(&tptr3,&nn3);CHKERRQ(info);

468:   for (i=0;i<nn1;i++) tptr1[i] = TaoMax( tptr2[i] , tptr3[i]);

470:   info = this->RestoreFloats(&tptr1,&nn1);CHKERRQ(info);
471:   info = vv->RestoreFloats(&tptr2,&nn2);CHKERRQ(info);
472:   info = ww->RestoreFloats(&tptr3,&nn3);CHKERRQ(info);

474:   TaoFunctionReturn(0);
475: }

479: int TaoVecFloatArray::Waxpby  ( double a, TaoVec* tx, double b, TaoVec* ty){
480:   TaoInt i,nn1,nn2,nn3;
481:   int info;
482:   float *tptr1,*tptr2,*tptr3;
483:   TaoVecFloatArray* xx =  (TaoVecFloatArray*)(tx);
484:   TaoVecFloatArray* yy =  (TaoVecFloatArray*)(ty);

486:   TaoFunctionBegin;
487:   info = this->GetFloats(&tptr1,&nn1);CHKERRQ(info);
488:   info = xx->GetFloats(&tptr2,&nn2);CHKERRQ(info);
489:   info = yy->GetFloats(&tptr3,&nn3);CHKERRQ(info);
490:   if (nn1!=nn2 || nn2!=nn3) {TaoFunctionReturn(1);}
491:   for (i=0;i<nn1;i++){ tptr1[i] = a * tptr2[i] + b * tptr3[i]; }
492:   info = this->RestoreFloats(&tptr1,&nn1);CHKERRQ(info);
493:   info = xx->RestoreFloats(&tptr2,&nn2);CHKERRQ(info);
494:   info = yy->RestoreFloats(&tptr3,&nn3);CHKERRQ(info);
495:   TaoFunctionReturn(0);
496: }

500: int TaoVecFloatArray::AbsoluteValue(){
501:   TaoInt i,nn;
502:   int info;
503:   float *vv;

505:   TaoFunctionBegin;
506:   info = this->GetFloats(&vv,&nn);CHKERRQ(info);
507:   for (i=0;i<nn;i++){ 
508:     vv[i]= TaoAbsScalar(vv[i]);
509:   }
510:   info = this->RestoreFloats(&vv,&nn);CHKERRQ(info);
511:   TaoFunctionReturn(0);
512: }

516: int TaoVecFloatArray::MinElement(double*val){
517:   TaoInt i,nn;
518:   int info;
519:   float dd=TAO_INFINITY,*vv;

521:   TaoFunctionBegin;
522:   info = this->GetFloats(&vv,&nn);CHKERRQ(info);
523:   for (i=0;i<nn;i++){ 
524:     if (vv[i]<dd) dd=vv[i];
525:   }
526:   info = this->RestoreFloats(&vv,&nn);CHKERRQ(info);
527:   *val = dd;
528:   TaoFunctionReturn(0);
529: }

533: int TaoVecFloatArray::StepMax(TaoVec*tv,double*step1){
534:   TaoInt i,nn1,nn2;
535:   int info;
536:   float *xx,*dx;
537:   float stepmax1=TAO_INFINITY;
538:   TaoVecFloatArray* vv =  (TaoVecFloatArray*)(tv);

540:   TaoFunctionBegin;
541:   info = this->GetFloats(&xx,&nn1);CHKERRQ(info);
542:   info = vv->GetFloats(&dx,&nn2);CHKERRQ(info);
543:   
544:   for (i=0;i<nn1;i++){
545:     if (xx[i] < 0){
546:       TaoFunctionReturn(1);
547:     } else if (dx[i]<0){ 
548:       stepmax1=TaoMin(stepmax1,-xx[i]/dx[i]);
549:     }
550:   }
551:   *step1=stepmax1;

553:   info = vv->RestoreFloats(&dx,&nn2);CHKERRQ(info);
554:   info = this->RestoreFloats(&xx,&nn1);CHKERRQ(info);
555:   TaoFunctionReturn(0);
556: }

560: int TaoVecFloatArray::StepMax2(TaoVec*tv,TaoVec*txl,TaoVec*txu, double*step2){
561:   TaoInt i,nn1,nn2;
562:   int info;
563:   float *xx,*dx,*xl,*xu;
564:   float stepmax2=0;
565:   TaoVecFloatArray* vv =  (TaoVecFloatArray*)(tv);
566:   TaoVecFloatArray* xxll =  (TaoVecFloatArray*)(txl);
567:   TaoVecFloatArray* xxuu =  (TaoVecFloatArray*)(txu);

569:   TaoFunctionBegin;
570:   info = this->GetFloats(&xx,&nn1);CHKERRQ(info);
571:   info = vv->GetFloats(&dx,&nn2);CHKERRQ(info);
572:   info = xxll->GetFloats(&xl,&nn2);CHKERRQ(info);
573:   info = xxuu->GetFloats(&xu,&nn2);CHKERRQ(info);
574:   
575:   for (i=0;i<nn1;i++){
576:     if (dx[i] > 0){
577:       stepmax2=TaoMax(stepmax2,(xu[i]-xx[i])/dx[i]);      
578:     } else if (dx[i]<0){ 
579:       stepmax2=TaoMax(stepmax2,(xl[i]-xx[i])/dx[i]);
580:     }
581:   }
582:   info = this->RestoreFloats(&xx,&nn1);CHKERRQ(info);
583:   info = vv->RestoreFloats(&dx,&nn2);CHKERRQ(info);
584:   info = xxll->RestoreFloats(&xl,&nn2);CHKERRQ(info);
585:   info = xxuu->RestoreFloats(&xu,&nn2);CHKERRQ(info);
586:   
587:   *step2=stepmax2;
588:   TaoFunctionReturn(0);
589: }


594: int TaoVecFloatArray::BoundGradientProjection(TaoVec*tg,TaoVec*txl,TaoVec*tx, TaoVec*txu){
595:   TaoInt i,nn1,nn2,nn3,nn4,nn5;
596:   int info;
597:   float *xptr,*xlptr,*xuptr,*gptr,*gpptr;
598:   TaoVecFloatArray* gg =  (TaoVecFloatArray*)(tg);
599:   TaoVecFloatArray* xxll =  (TaoVecFloatArray*)(txl);
600:   TaoVecFloatArray* xx =  (TaoVecFloatArray*)(tx);
601:   TaoVecFloatArray* xxuu =  (TaoVecFloatArray*)(txu);

603:   TaoFunctionBegin;
604:   info = this->GetFloats(&gpptr,&nn1);CHKERRQ(info);
605:   info = gg->GetFloats(&gptr,&nn2);CHKERRQ(info);
606:   info = xxll->GetFloats(&xlptr,&nn3);CHKERRQ(info);
607:   info = xx->GetFloats(&xptr,&nn4);CHKERRQ(info);
608:   info = xxuu->GetFloats(&xuptr,&nn5);CHKERRQ(info);
609:   if (nn1!=nn2 || nn2!=nn3 || nn3!=nn4 || nn4!=nn5) {TaoFunctionReturn(1);}

611:   for (i=0; i<nn1; i++){

613:     gpptr[i] = gptr[i];
614:     if (gpptr[i]>0 && xptr[i]<=xlptr[i]){
615:       gpptr[i] = 0;
616:     } else if (gpptr[i]<0 && xptr[i]>=xuptr[i]){
617:       gpptr[i] = 0;
618:     }
619:   }
620:   info = this->RestoreFloats(&gpptr,&nn1);CHKERRQ(info);
621:   info = gg->RestoreFloats(&gptr,&nn2);CHKERRQ(info);
622:   info = xxll->RestoreFloats(&xlptr,&nn3);CHKERRQ(info);
623:   info = xx->RestoreFloats(&xptr,&nn4);CHKERRQ(info);
624:   info = xxuu->RestoreFloats(&xuptr,&nn5);CHKERRQ(info);

626:   TaoFunctionReturn(0);
627: }

629: inline static float fischer(float a, float b) 
630: {

632: #ifdef TODD

634:   if (TaoAbsScalar(a) > TaoAbsScalar(b)) {
635:     return sqrt(a*a + b*b) - a - b;
636:   }
637:   return sqrt(a*a + b*b) - b - a;

639: #else

641:    // Method suggested by Bob Vanderbei

643:    if (a + b <= 0) {
644:      return sqrt(a*a + b*b) - (a + b);
645:    }
646:    return -2.0*a*b / (sqrt(a*a + b*b) + (a + b));

648: #endif

650: }

654: int TaoVecFloatArray::Fischer(TaoVec *tx, TaoVec *tf, TaoVec *tl, TaoVec *tu)
655: {
656:   TaoInt i,nn1,nn2,nn3,nn4,nn5;
657:   int info;
658:   float *vv,*x,*f,*l,*u;

660:   TaoVecFloatArray* xx =  (TaoVecFloatArray*)(tx);
661:   TaoVecFloatArray* ff =  (TaoVecFloatArray*)(tf);
662:   TaoVecFloatArray* ll =  (TaoVecFloatArray*)(tl);
663:   TaoVecFloatArray* uu =  (TaoVecFloatArray*)(tu);

665:   TaoFunctionBegin;
666:   info = this->GetFloats(&vv,&nn1);CHKERRQ(info);
667:   info = xx->GetFloats(&x,&nn2);CHKERRQ(info);
668:   info = ff->GetFloats(&f,&nn3);CHKERRQ(info);
669:   info = ll->GetFloats(&l,&nn4);CHKERRQ(info);
670:   info = uu->GetFloats(&u,&nn5);CHKERRQ(info);
671:   if (nn1!=nn2 || nn2!=nn3 || nn3!=nn4 || nn4!=nn5) {
672:     TaoFunctionReturn(1);
673:   }

675:   for (i=0;i<nn1;i++) {

677:     if ((l[i] <= -TAO_INFINITY) && (u[i] >= TAO_INFINITY)) {
678:       vv[i] = -f[i];
679:     } 
680:     else if (l[i] <= -TAO_INFINITY) {
681:       vv[i] = -fischer(u[i] - x[i], -f[i]);
682:     } 
683:     else if (u[i] >=  TAO_INFINITY) {
684:       vv[i] =  fischer(x[i] - l[i],  f[i]);
685:     } 
686:     else if (l[i] == u[i]) {
687:       vv[i] = l[i] - x[i];
688:     }
689:     else {
690:       vv[i] =  fischer(u[i] - x[i], -f[i]);
691:       vv[i] =  fischer(x[i] - l[i],  vv[i]);
692:     }

694:   }

696:   info = this->RestoreFloats(&vv,&nn1);CHKERRQ(info);
697:   info = xx->RestoreFloats(&x,&nn2);CHKERRQ(info);
698:   info = ff->RestoreFloats(&f,&nn3);CHKERRQ(info);
699:   info = ll->RestoreFloats(&l,&nn4);CHKERRQ(info);
700:   info = uu->RestoreFloats(&u,&nn5);CHKERRQ(info);

702:   TaoFunctionReturn(0);
703: }

705: inline static float sfischer(float a, float b, float c)
706: {

708: #ifdef TODD

710:   if (TaoAbsScalar(a) > TaoAbsScalar(b)) {
711:     return sqrt(a*a + b*b + 2.0*c*c) - a - b;
712:   }
713:   return sqrt(a*a + b*b + 2.0*c*c) - b - a;

715: #else

717:    // Method suggested by Bob Vanderbei

719:    if (a + b <= 0) {
720:      return sqrt(a*a + b*b + 2.0*c*c) - (a + b);
721:    }
722:    return 2.0*(c*c - a*b) / (sqrt(a*a + b*b + 2.0*c*c) + (a + b));

724: #endif

726: }

730: int TaoVecFloatArray::SFischer(TaoVec *tx, TaoVec *tf, TaoVec *tl, TaoVec *tu, double mu)
731: {
732:   TaoInt i, nn1, nn2, nn3, nn4, nn5;
733:   int info;
734:   float *vv, *x, *f, *l, *u;

736:   TaoVecFloatArray *xx = (TaoVecFloatArray *)(tx);
737:   TaoVecFloatArray *ff = (TaoVecFloatArray *)(tf);
738:   TaoVecFloatArray *ll = (TaoVecFloatArray *)(tl);
739:   TaoVecFloatArray *uu = (TaoVecFloatArray *)(tu);

741:   TaoFunctionBegin;

743:   if ((mu >= -TAO_EPSILON) && (mu <= TAO_EPSILON)) {
744:     Fischer(tx, tf, tl, tu);
745:   }
746:   else {
747:     info = this->GetFloats(&vv, &nn1); CHKERRQ(info);
748:     info = xx->GetFloats(&x, &nn2); CHKERRQ(info);
749:     info = ff->GetFloats(&f, &nn3); CHKERRQ(info);
750:     info = ll->GetFloats(&l, &nn4); CHKERRQ(info);
751:     info = uu->GetFloats(&u, &nn5); CHKERRQ(info);

753:     if (nn1!=nn2 || nn2!=nn3 || nn3!=nn4 || nn4!=nn5) {
754:       TaoFunctionReturn(1);
755:     }

757:     for (i = 0; i < nn1; ++i) {
758:       if ((l[i] <= -TAO_INFINITY) && (u[i] >= TAO_INFINITY)) {
759:         vv[i] = -f[i] - mu*x[i];
760:       }
761:       else if (l[i] <= -TAO_INFINITY) {
762:         vv[i] = -sfischer(u[i] - x[i], -f[i], mu);
763:       }
764:       else if (u[i] >=  TAO_INFINITY) {
765:         vv[i] =  sfischer(x[i] - l[i],  f[i], mu);
766:       }
767:       else if (l[i] == u[i]) {
768:         vv[i] = l[i] - x[i];
769:       }
770:       else {
771:         vv[i] =  sfischer(u[i] - x[i], -f[i], mu);
772:         vv[i] =  sfischer(x[i] - l[i],  vv[i], mu);
773:       }
774:     }

776:     info = this->RestoreFloats(&vv, &nn1); CHKERRQ(info);
777:     info = xx->RestoreFloats(&x, &nn2); CHKERRQ(info);
778:     info = ff->RestoreFloats(&f, &nn3); CHKERRQ(info);
779:     info = ll->RestoreFloats(&l, &nn4); CHKERRQ(info);
780:     info = uu->RestoreFloats(&u, &nn5); CHKERRQ(info);
781:   }
782:   TaoFunctionReturn(0);
783: }