Actual source code: taovec.c

  1: #include "tao_general.h"          /*I "tao_solver.h"  I*/
  2: #include "taovec.h"

  6: /*@C
  7:    TaoVecDestroy - Destroys the TaoVec object.

  9:    Input Parameter:
 10: .  vv - the vector

 12:    Level: beginner
 13: @*/
 14: int TaoVecDestroy( TaoVec* vv){
 15:   TaoFunctionBegin;
 16:   if (vv!=TAO_NULL && vv!=0){ 
 17:     delete vv;
 18:   }
 19:   vv=TAO_NULL;
 20:   TaoFunctionReturn(0);
 21: }

 25: /*@C
 26:    CloneVecs - Creates an array of pointers to new TaoVec objects. The new 
 27:    objects have the same structure as this one.

 29:    Input Parameter:
 30: .  nn -  number of new vectors

 32:    Output Parameter:
 33: .  tvs -  pointer to array TaoVec pointers

 35: .seealso TaoVec::Clone()

 37:    Level: intermediate
 38: @*/
 39: int TaoVec::CloneVecs(TaoInt nn, TaoVec***tvs){
 40:   int info;
 41:   TaoInt i;
 42:   TaoVec ** ntv;
 43:   TaoFunctionBegin;
 44:   ntv = new TaoVec*[nn];
 45:   for (i=0;i<nn;i++){
 46:     info = this->Clone(&ntv[i]);CHKERRQ(info);
 47:   }
 48:   *tvs=ntv;
 49:   TaoFunctionReturn(0);
 50: }

 54: /*@C
 55:    DestroyVecs - Destroys an array TaoVec objects.

 57:    Input Parameter:
 58: +  nn -  number of new vectors
 59: -  tvs -  pointer to array TaoVec pointers

 61:    Level: advanced

 63: .seealso TaoVec::CloneVecs()
 64: @*/
 65: int TaoVec::DestroyVecs(TaoInt nn, TaoVec**ntv){
 66:   int info;
 67:   TaoInt i;
 68:   TaoFunctionBegin;
 69:   for (i=0;i<nn;i++){
 70:     info = TaoVecDestroy( ntv[i] );CHKERRQ(info);
 71:   }
 72:   delete[] ntv;
 73:   TaoFunctionReturn(0);
 74: }

 78: /*@C
 79:    Clone - Creates a new TaoVec object with the same structure as this one.  It does not
 80:    copy the value to the new vector.

 82:    Input:
 83: .  vv - address of a pointer to a TaoVec

 85:    Output Parameter:
 86: .  vv - address of a pointer to new TaoVec object

 88: .seealso TaoVec::CloneVecs(), TaoVec::CopyFrom()

 90:    Level: intermediate
 91: @*/
 92: int TaoVec::Clone( TaoVec* *vv ){
 93:   TaoFunctionBegin;
 94:   SETERRQ(56,"Operation not defined");
 95:   /* TaoFunctionReturn(1); */
 96: }

100: /*@C
101:    Compatible - Determines whether this vector belongs to the same space as another,
102:    and operations such as inner product and sum are well defined.

104:    Input Parameter:
105: .  vv -  TAO vector to which to the comparison is made

107:    Output Value:
108: .  flag - TAO_TRUE if the two vectors are Compatible and TAO_FALSE otherwise.

110:    Level: advanced

112: .seealso TaoVec::GetDimension()
113: @*/
114: int TaoVec::Compatible(TaoVec* vv, TaoTruth *flag){
115:   TaoFunctionBegin;
116:   if (!flag){
117:     TaoFunctionReturn(1);
118:   }
119:   *flag=TAO_FALSE;
120:   TaoFunctionReturn(0);
121: }

125: /*@C
126:    SetToConstant - Sets each element of this vector equal to the specified constant.

128:    Input Parameter:
129: .  c -  a constant

131:    Level: intermediate

133: .seealso TaoVec::Scale()
134: @*/
135: int TaoVec::SetToConstant( double c ){
136:   int info;
137:   TaoInt nn,i;
138:   TaoScalar *tptr;

140:   TaoFunctionBegin;
141:   info = this->GetArray(&tptr,&nn);CHKERRQ(info);
142:   for (i=0;i<nn;i++){ tptr[i]=c; }
143:   info = this->RestoreArray(&tptr,&nn);CHKERRQ(info);
144:   TaoFunctionReturn(0);
145: }

149: /*@C
150:    SetToZero - Sets each element of this vector equal to zero.

152:    Input Parameters: none

154:    Level: intermediate

156: .seealso TaoVec::SetToConstant()
157: @*/
158: int TaoVec::SetToZero(){
159:   int info;
160:   TaoInt nn,i;
161:   TaoScalar *tptr;

163:   TaoFunctionBegin;
164:   info = this->GetArray(&tptr,&nn);CHKERRQ(info);
165:   for (i=0;i<nn;i++){ tptr[i]=0; }
166:   info = this->RestoreArray(&tptr,&nn);CHKERRQ(info);
167:   TaoFunctionReturn(0);
168: }

172: /*@C
173:    CopyFrom - Copies the contents of one vector into this vector.

175:    Input Parameter:
176: .  vv -  A TaoVec from which the contents will be copied.

178:    Level: intermediate

180: .seealso TaoVec::Axpy(), TaoVec::ScaleCopyFrom()
181: @*/
182: int TaoVec::CopyFrom( TaoVec* vv ){
183:   int info;
184:   TaoInt nn1,nn2,i;
185:   TaoScalar *tptr1,*tptr2;

187:   TaoFunctionBegin;
188:   info = this->GetArray(&tptr1,&nn1);CHKERRQ(info);
189:   if (vv!=this){
190:     info = vv->GetArray(&tptr2,&nn2);CHKERRQ(info);
191:     if (nn1!=nn2) {TaoFunctionReturn(1);}
192:     for (i=0;i<nn1;i++){ tptr1[i]=tptr2[i]; }
193:     info = vv->RestoreArray(&tptr2,&nn2);CHKERRQ(info);
194:   }
195:   info = this->RestoreArray(&tptr1,&nn1);CHKERRQ(info);
196:   TaoFunctionReturn(0);
197: }

201: /*@C
202:    ScaleCopyFrom - Copies the contents of one vector into this vector and scales it.

204:    Input Parameter:
205: +  a - the scalar 
206: -  vv -  A TaoVec from which the contents will be copied.

208:    Level: intermediate

210: .seealso TaoVec::Axpy(), TaoVec::Aypx()
211: @*/
212: int TaoVec::ScaleCopyFrom( double a, TaoVec* vv ){
213:   int info;
214:   TaoInt nn1,nn2,i;
215:   TaoScalar *tptr1,*tptr2;

217:   TaoFunctionBegin;
218:   info = this->GetArray(&tptr1,&nn1);CHKERRQ(info);
219:   if (vv!=this){
220:     info = vv->GetArray(&tptr2,&nn2);CHKERRQ(info);
221:     if (nn1!=nn2) {TaoFunctionReturn(1);}
222:     for (i=0;i<nn1;i++){ tptr1[i]=a*tptr2[i]; }
223:     info = vv->RestoreArray(&tptr2,&nn2);CHKERRQ(info);
224:   } else {
225:     this->Scale(a);
226:   }
227:   info = this->RestoreArray(&tptr1,&nn1);CHKERRQ(info);
228:   TaoFunctionReturn(0);
229: }

233: /*@C
234:    NormInfinity - Computes the infinity norm of this vector.

236:    Ouput Parameter:
237: .  vnorm -  the infinity norm of this vector

239:    Level: intermediate

241: .seealso TaoVec::Norm1(), TaoVec::Norm2()
242: @*/
243: int TaoVec::NormInfinity(double *vnorm){
244:   int  info;
245:   TaoInt nn,i;
246:   TaoScalar dd=0, *v;

248:   TaoFunctionBegin;
249:   info = this->GetArray(&v,&nn);CHKERRQ(info);
250:   for (i=0;i<nn;i++){
251:     if (v[i]<0 && v[i]<-dd) dd=-v[i];
252:     else if (v[i]>0 && v[i]>dd) dd=v[i];
253:   }
254:   info = this->RestoreArray(&v,&nn);CHKERRQ(info);
255:   *vnorm=dd;
256:   TaoFunctionReturn(0);
257: }

261: /*@C
262:    Norm1 - Computes the one-norm of this vector.

264:    Ouput Parameter:
265: .  vnorm -  the one-norm of this vector

267:    Level: intermediate

269: .seealso TaoVec::NormInfinity(), TaoVec::Norm2()
270: @*/
271: int TaoVec::Norm1(double *vnorm){
272:   int info;
273:   TaoInt nn,i;
274:   TaoScalar dd=0, *v;

276:   TaoFunctionBegin;
277:   info = this->GetArray(&v,&nn);CHKERRQ(info);
278:   for (i=0;i<nn;i++){ 
279:     if (v[i]<0) dd-=v[i]; 
280:     else dd+=v[i];
281:   }
282:   info = this->RestoreArray(&v,&nn);CHKERRQ(info);
283:   *vnorm = dd;
284:   TaoFunctionReturn(0);
285: }

289: /*@C
290:    Norm2 - Compute the two-norm of this vector.

292:    Ouput Parameter:
293: .  vnorm -  the two-norm of this vector

295:    Level: intermediate

297: .seealso TaoVec::Norm1(), TaoVec::NormInfinity(), TaoVec::Norm2squared()
298: @*/
299: int TaoVec::Norm2(double *vnorm){
300:   int info;
301:   TaoInt nn,i;
302:   TaoScalar dd=0, *v;

304:   TaoFunctionBegin;
305:   info = this->GetArray(&v,&nn);CHKERRQ(info);
306:   for (i=0;i<nn;i++) dd+=v[i]*v[i];
307:   *vnorm=sqrt(dd);
308:   info = this->RestoreArray(&v,&nn);CHKERRQ(info);
309:   TaoFunctionReturn(0);
310: }

314: /*@C
315:    Norm2squared - Computes the square of the two norm of this vector.

317:    Ouput Parameter:
318: .  vnorm2 -  the square of the two norm of this vector

320:    Level: intermediate

322: .seealso TaoVec::Norm2()
323: @*/
324: int TaoVec::Norm2squared(double *vnorm2){
325:   int info;
326:   TaoInt nn,i;
327:   TaoScalar dd=0, *v;

329:   TaoFunctionBegin;
330:   info = this->GetArray(&v,&nn);CHKERRQ(info);
331:   for (i=0;i<nn;i++) dd+=v[i]*v[i];
332:   *vnorm2=dd;
333:   info = this->RestoreArray(&v,&nn);CHKERRQ(info);
334:   TaoFunctionReturn(0);
335: }

339: /*@C
340:    Scale - Multiplies this vector by a scalar.

342:    Input Parameter:
343: .  alpha -  the scalar

345:    Level: intermediate

347: .seealso TaoVec::SetToConstant(),  TaoVec::Aypx()
348: @*/
349: int TaoVec::Scale( double alpha ){
350:   int info;
351:   TaoInt nn,i;
352:   TaoScalar *v;

354:   TaoFunctionBegin;
355:   info = this->GetArray(&v,&nn);CHKERRQ(info);
356:   for (i=0;i<nn;i++) v[i]*=alpha;
357:   info = this->RestoreArray(&v,&nn);CHKERRQ(info);
358:   TaoFunctionReturn(0);
359: }

363: /*@C
364:    Axpy - Adds a scalar multiple of a vector to this vector. (this += alpha * xx)

366:    Input Parameter:
367: +  alpha -  the scalar 
368: -  xx - the vector

370:    Level: intermediate

372: .seealso TaoVec::CopyFrom(), TaoVec::Aypx() 
373: @*/
374: int TaoVec::Axpy( double alpha, TaoVec* xx ){
375:   int info;
376:   TaoInt nn1,nn2,i;
377:   TaoScalar *tptr1,*tptr2;

379:   TaoFunctionBegin;
380:   info = this->GetArray(&tptr1,&nn1);CHKERRQ(info);
381:   info = xx->GetArray(&tptr2,&nn2);CHKERRQ(info);
382:   if (nn1!=nn2) {TaoFunctionReturn(1);}
383:   for (i=0;i<nn1;i++){ tptr1[i]+= alpha * tptr2[i]; }
384:   info = this->RestoreArray(&tptr1,&nn1);CHKERRQ(info);
385:   info = xx->RestoreArray(&tptr2,&nn2);CHKERRQ(info);
386:   TaoFunctionReturn(0);
387: }

391: /*@C
392:    Axpby - Adds a scalar multiple of a vector to a multiple of this vector. (this=alpha*xx + beta*this)

394:    Input Parameter:
395: +  alpha -  the scalar of tx
396: .  xx - the vector
397: -  beta -  the scalar multiple of this vector

399:    Level: intermediate

401: .seealso TaoVec::Axpy(), TaoVec::Aypx() 
402: @*/
403: int TaoVec::Axpby( double alpha, TaoVec* xx, double beta ){
404:   int info;
405:   TaoInt nn1,nn2,i;
406:   TaoScalar *tptr1,*tptr2;

408:   TaoFunctionBegin;
409:   info = this->GetArray(&tptr1,&nn1);CHKERRQ(info);
410:   info = xx->GetArray(&tptr2,&nn2);CHKERRQ(info);
411:   if (nn1!=nn2) {TaoFunctionReturn(1);}
412:   for (i=0;i<nn1;i++){ tptr1[i] = beta * tptr1[i] + alpha * tptr2[i]; }
413:   info = this->RestoreArray(&tptr1,&nn1);CHKERRQ(info);
414:   info = xx->RestoreArray(&tptr2,&nn2);CHKERRQ(info);
415:   TaoFunctionReturn(0);
416: }

420: /*@C
421:    Aypx - Adds a vector to a scalar multiple of this vector. (this=alpha*this+xx)

423:    Input Parameter:
424: +  alpha -  the scalar 
425: -  xx - the vector

427:    Level: intermediate

429: .seealso TaoVec::Scale(), TaoVec::Axpy()
430: @*/
431: int TaoVec::Aypx( double alpha, TaoVec* xx ){
432:   int info;
433:   TaoInt nn1,nn2,i;
434:   TaoScalar *tptr1,*tptr2;

436:   TaoFunctionBegin;
437:   info = this->GetArray(&tptr1,&nn1);CHKERRQ(info);
438:   info = xx->GetArray(&tptr2,&nn2);CHKERRQ(info);
439:   if (nn1!=nn2) {TaoFunctionReturn(1);}
440:   for (i=0;i<nn1;i++){ tptr1[i] = alpha * tptr1[i] + tptr2[i]; }
441:   info = this->RestoreArray(&tptr1,&nn1);CHKERRQ(info);
442:   info = xx->RestoreArray(&tptr2,&nn2);CHKERRQ(info);
443:   TaoFunctionReturn(0);
444: }
445:  
448: /*@C
449:    AddConstant - Adds a constant to each element of this vector.

451:    Input Parameter:
452: .  alpha -  the scalar 

454:    Level: intermediate

456: .seealso TaoVec::SetToConstant(), TaoVec::Axpy() 
457: @*/
458: int TaoVec::AddConstant( double alpha ){
459:   int info;
460:   TaoInt nn,i;
461:   TaoScalar *v;

463:   TaoFunctionBegin;
464:   info = this->GetArray(&v,&nn);CHKERRQ(info);
465:   for (i=0;i<nn;i++) v[i]+=alpha;
466:   info = this->RestoreArray(&v,&nn);CHKERRQ(info);
467:   TaoFunctionReturn(0);
468: }

472: /*@C
473:    Dot - Computes the inner product of this vector with another vector.

475:    Input Parameter:
476: .  vv -  another TaoVec object

478:    Output Parameter:
479: .  vDotv -  the inner product of the two vectors

481:    Level: intermediate

483: .seealso TaoVec::Norm() 
484: @*/
485: int TaoVec::Dot( TaoVec* vv, double *vDotv ){
486:   int info;
487:   TaoInt nn1,nn2,i;
488:   TaoScalar dd=0,*tptr1,*tptr2;

490:   TaoFunctionBegin;
491:   info = this->GetArray(&tptr1,&nn1);CHKERRQ(info);
492:   info = vv->GetArray(&tptr2,&nn2);CHKERRQ(info);
493:   if (nn1!=nn2) {TaoFunctionReturn(1);}
494:   for (i=0;i<nn1;i++) dd+=tptr1[i]*tptr2[i];
495:   info = this->RestoreArray(&tptr1,&nn1);CHKERRQ(info);
496:   info = vv->RestoreArray(&tptr2,&nn2);CHKERRQ(info);
497:   *vDotv=dd;
498:   TaoFunctionReturn(0);
499: }

503: /*@C
504:    Negate - Multiplies the elements of this vector by negative one.

506:    Input Parameters: none

508:    Level: intermediate

510: .seealso TaoVec::Scale()
511: @*/
512: int TaoVec::Negate(){ 
513:   int info;
514:   TaoInt nn,i;
515:   TaoScalar *v;

517:   TaoFunctionBegin;
518:   info = this->GetArray(&v,&nn);CHKERRQ(info);
519:   for (i=0;i<nn;i++) v[i]=-v[i];
520:   info = this->RestoreArray(&v,&nn);CHKERRQ(info);
521:   TaoFunctionReturn(0);
522: }

526: /*@C
527:    Reciprocal - Sets each element of this vector to its Reciprocal.

529:    Input Parameters: none

531:    Level: intermediate

533: .seealso TaoVec::PointwiseDivide()
534: @*/
535: int TaoVec::Reciprocal(){ 
536:   int info;
537:   TaoInt nn,i;
538:   TaoScalar *v;

540:   TaoFunctionBegin;
541:   info = this->GetArray(&v,&nn);CHKERRQ(info);
542:   for (i=0;i<nn;i++){ 
543:     if (v[i]!=0)  v[i]= 1.0/v[i];
544:     else v[i]=TAO_INFINITY;
545:   }
546:   info = this->RestoreArray(&v,&nn);CHKERRQ(info);
547:   TaoFunctionReturn(0);
548: }

552: /*@C
553:    Sqrt - Sets each element of this vector to its square root.

555:    Input Parameters: none

557:    Level: intermediate

559: @*/
560: int TaoVec::Sqrt(){ 
561:   int info;
562:   TaoInt i;
563:   TaoInt nn;
564:   TaoScalar *v;

566:   TaoFunctionBegin;
567:   info = this->GetArray(&v,&nn);CHKERRQ(info);
568:   for (i=0;i<nn;i++){ 
569:     if (v[i] >= 0) {
570:       v[i] = sqrt(v[i]);
571:     }
572:     else {
573:       v[i] = TAO_INFINITY;
574:     }
575:   }
576:   info = this->RestoreArray(&v,&nn);CHKERRQ(info);
577:   TaoFunctionReturn(0);
578: }

582: /*@C
583:    Pow - Raises each element of this vector to a power.

585:    Input Parameters: none

587:    Level: intermediate

589: @*/
590: int TaoVec::Pow(double p){ 
591:   int info;
592:   TaoInt nn,i;
593:   TaoScalar *v;

595:   TaoFunctionBegin;
596:   info = this->GetArray(&v,&nn);CHKERRQ(info);
597:   for (i=0;i<nn;i++){ 
598:     if (v[i] >= 0) {
599:       v[i] = pow(v[i], p);
600:     }
601:     else {
602:       v[i] = TAO_INFINITY;
603:     }
604:   }
605:   info = this->RestoreArray(&v,&nn);CHKERRQ(info);
606:   TaoFunctionReturn(0);
607: }

611: /*@C
612:    GetDimension - Gets the dimension of the vector space where this vector belongs.

614:    Output Parameter:
615: .  n - the dimension of the vector space

617:    Level: intermediate

619: .seealso TaoVec::Compatible()
620: @*/
621: int TaoVec::GetDimension(TaoInt *n){
622:   TaoFunctionBegin;
623:   SETERRQ(56,"Operation not defined");
624:   /* TaoFunctionReturn(1); */
625: }

629: /*@C
630:    PointwiseMultiply - Computes the componentwise multiplication of two vectors 
631:    and stores the result in this vector.

633:    Input Parameters:
634: .   vv, ww - the two vectors

636:    Level: intermediate

638: .seealso TaoVec::PointwiseDivide()
639: @*/
640: int TaoVec::PointwiseMultiply( TaoVec* vv, TaoVec* ww ){
641:   int info;
642:   TaoInt nn1,nn2,nn3,i;
643:   TaoScalar *tptr1,*tptr2,*tptr3;

645:   TaoFunctionBegin;
646:   info = this->GetArray(&tptr1,&nn1);CHKERRQ(info);
647:   if (this!=vv){ 
648:     info = vv->GetArray(&tptr2,&nn2);CHKERRQ(info);
649:   } else {
650:     tptr2=tptr1; nn2=nn1;
651:   }
652:   if (this!=ww && vv!=ww){ 
653:     info = ww->GetArray(&tptr3,&nn3);CHKERRQ(info);
654:   } else if (vv==ww){ 
655:     tptr3=tptr2; nn3=nn2;
656:   } else {
657:     tptr3=tptr1; nn3=nn1;
658:   }
659:   if (nn1!=nn2 || nn2!=nn3) {TaoFunctionReturn(1);}
660:   for (i=0;i<nn1;i++) tptr1[i]=tptr2[i] * tptr3[i];
661:   info = this->RestoreArray(&tptr1,&nn1);CHKERRQ(info);
662:   if (this!=vv){ 
663:     info = vv->RestoreArray(&tptr2,&nn2);CHKERRQ(info);
664:   }
665:   if (this!=ww && vv!=ww){ 
666:     info = ww->RestoreArray(&tptr3,&nn3);CHKERRQ(info);
667:   }

669:   TaoFunctionReturn(0);
670: }

674: /*@C
675:    PointwiseDivide - Computes the componentwise division of two vectors 
676:    and stores the result in this vector.

678:    Input Parameters:
679: +  vv - the vector of numerators
680: -  ww - the vector of denominators

682:    Level: intermediate

684: .seealso TaoVec::PointwiseMultiply()
685: @*/
686: int TaoVec::PointwiseDivide( TaoVec* vv , TaoVec* ww){
687:   int info;
688:   TaoInt nn1,nn2,nn3,i;
689:   TaoScalar *tptr1,*tptr2,*tptr3;

691:   TaoFunctionBegin;
692:   info = this->GetArray(&tptr1,&nn1);CHKERRQ(info);
693:   if (this!=vv){ 
694:     info = vv->GetArray(&tptr2,&nn2);CHKERRQ(info);
695:   } else {
696:     tptr2=tptr1; nn2=nn1;
697:   }
698:   if (this!=ww && vv!=ww){ 
699:     info = ww->GetArray(&tptr3,&nn3);CHKERRQ(info);
700:   } else if (vv==ww){
701:     tptr3=tptr2; nn3=nn2;
702:   } else {
703:     tptr3=tptr1; nn3=nn1;
704:   }
705:   if (nn1!=nn2 || nn2!=nn3) {TaoFunctionReturn(1);}
706:   for (i=0;i<nn1;i++) tptr1[i]=tptr2[i] / tptr3[i];
707:   info = this->RestoreArray(&tptr1,&nn1);CHKERRQ(info);
708:   if (this!=vv){ 
709:     info = vv->RestoreArray(&tptr2,&nn2);CHKERRQ(info);
710:   }
711:   if (this!=ww && vv!=ww){ 
712:     info = ww->RestoreArray(&tptr3,&nn3);CHKERRQ(info);
713:   }
714:   TaoFunctionReturn(0);
715: }

719: /*@C
720:    Median - Computes the componentwise median of three vectors 
721:    and stores the result in this vector.

723:    Input Parameters:
724: .   vv, ww, xx - the three vectors

726:    Level: intermediate

728: @*/
729: int TaoVec::Median( TaoVec* vv, TaoVec* ww, TaoVec* xx){
730:   int info;
731:   TaoInt nn1,nn2,nn3,nn4,i;
732:   TaoScalar *tptr1,*tptr2,*tptr3,*tptr4;

734:   TaoFunctionBegin;
735:   info = this->GetArray(&tptr1,&nn1);CHKERRQ(info);
736:   if (this!=vv){
737:     info = vv->GetArray(&tptr2,&nn2);CHKERRQ(info);
738:   } else {
739:     tptr2=tptr1; nn2=nn1;
740:   }
741:   if (this!=ww){
742:     info = ww->GetArray(&tptr3,&nn3);CHKERRQ(info);
743:   } else {
744:     tptr3=tptr1; nn3=nn1;
745:   }
746:   if (this!=xx){
747:     info = xx->GetArray(&tptr4,&nn4);CHKERRQ(info);
748:   } else {
749:     tptr4=tptr1; nn4=nn1;
750:   }

752:   if (nn1!=nn2 || nn2!=nn3 || nn3!=nn4) {TaoFunctionReturn(1);}
753:   for (i=0;i<nn1;i++){
754:     if (tptr2[i]<=tptr3[i] && tptr3[i] <= tptr4[i]){
755:       tptr1[i]=tptr3[i];
756:     } else if (tptr4[i]<=tptr3[i] && tptr3[i] <= tptr2[i]){
757:       tptr1[i]=tptr3[i];
758:     } else if (tptr3[i]<=tptr2[i] && tptr2[i] <= tptr4[i]){
759:       tptr1[i]=tptr2[i];
760:     } else if (tptr4[i]<=tptr2[i] && tptr2[i] <= tptr3[i]){
761:       tptr1[i]=tptr2[i];
762:     } else {
763:       tptr1[i]=tptr4[i];
764:     }
765:   }
766:   info = this->RestoreArray(&tptr1,&nn1);CHKERRQ(info);
767:   if (this!=vv){
768:     info = vv->RestoreArray(&tptr2,&nn2);CHKERRQ(info);
769:   }
770:   if (this!=ww){
771:     info = ww->RestoreArray(&tptr3,&nn3);CHKERRQ(info);
772:   }
773:   if (this!=xx){
774:     info = xx->RestoreArray(&tptr4,&nn4);CHKERRQ(info);
775:   }
776:   TaoFunctionReturn(0);
777: }

781: /*@C
782:    PointwiseMinimum - Computes the componentwise minimum of two vectors 
783:    and stores the result in this vector.

785:    Input Parameters:
786: .   vv, ww - the two vectors

788:    Level: intermediate

790: .seealso TaoVec::PointwiseMaximum()
791: @*/
792: int TaoVec::PointwiseMinimum( TaoVec* vv, TaoVec* ww){
793:   int info;
794:   TaoInt nn1,nn2,nn3,i;
795:   TaoScalar *tptr1,*tptr2,*tptr3;

797:   TaoFunctionBegin;
798:   info = this->GetArray(&tptr1,&nn1);CHKERRQ(info);
799:   if (vv!=this){
800:     info = vv->GetArray(&tptr2,&nn2);CHKERRQ(info);
801:   } else {
802:     tptr2=tptr1;  nn2=nn1;
803:   }
804:   if (ww!=this){
805:     info = ww->GetArray(&tptr3,&nn3);CHKERRQ(info);
806:   } else {
807:     tptr3=tptr1;  nn3=nn1;
808:   }

810:   if (nn1!=nn2 || nn2!=nn3) {TaoFunctionReturn(1);}
811:   for (i=0;i<nn1;i++) tptr1[i] = TaoMin( tptr2[i] , tptr3[i]);
812:   info = this->RestoreArray(&tptr1,&nn1);CHKERRQ(info);
813:   if (vv!=this){
814:     info = vv->RestoreArray(&tptr2,&nn2);CHKERRQ(info);
815:   }
816:   if (ww!=this){
817:     info = ww->RestoreArray(&tptr3,&nn3);CHKERRQ(info);
818:   }
819:   TaoFunctionReturn(0);
820: }

824: /*@C
825:    PointwiseMaximum - Computes the componentwise minimum of two vectors 
826:    and stores the result in this vector.

828:    Input Parameters:
829: .  vv, ww - the two vectors

831:    Level: intermediate

833: .seealso TaoVec::PointwiseMinimum()
834: @*/
835: int TaoVec::PointwiseMaximum( TaoVec* vv, TaoVec* ww){
836:   int info;
837:   TaoInt nn1,nn2,nn3,i;
838:   TaoScalar *tptr1,*tptr2,*tptr3;

840:   TaoFunctionBegin;
841:   info = this->GetArray(&tptr1,&nn1);CHKERRQ(info);
842:   if (vv!=this){
843:     info = vv->GetArray(&tptr2,&nn2);CHKERRQ(info);
844:   } else {
845:     tptr2=tptr1;  nn2=nn1;
846:   }
847:   if (ww!=this){
848:     info = ww->GetArray(&tptr3,&nn3);CHKERRQ(info);
849:   } else {
850:     tptr3=tptr1;  nn3=nn1;
851:   }
852:   if (nn1!=nn2 || nn2!=nn3) {TaoFunctionReturn(1);}
853:   for (i=0;i<nn1;i++) tptr1[i] = TaoMax( tptr2[i] , tptr3[i]);
854:   info = this->RestoreArray(&tptr1,&nn1);CHKERRQ(info);
855:   if (vv!=this){
856:     info = vv->RestoreArray(&tptr2,&nn2);CHKERRQ(info);
857:   }
858:   if (ww!=this){
859:     info = ww->RestoreArray(&tptr3,&nn3);CHKERRQ(info);
860:   }
861:   TaoFunctionReturn(0);
862: }

866: /*@C
867:    Waxpby - Sums two scaled vectors and stores the result in this vector. (this=alpha*xx+beta*yy)

869:    Input Parameters:
870: +  a - the multiple of the first vector
871: .  xx - the first vector
872: .  b - the multiple of the second vector
873: -  yy - the second vector

875:    Level: intermediate

877: .seealso TaoVec::Axpy(), TaoVec::Axpby(), TaoVec::Aypx();
878: @*/
879: int TaoVec::Waxpby  ( double a, TaoVec* xx, double b, TaoVec* yy){
880:   int info;
881:   TaoInt nn1,nn2,nn3,i;
882:   TaoScalar *tptr1,*tptr2,*tptr3;

884:   TaoFunctionBegin;
885:   info = this->GetArray(&tptr1,&nn1);CHKERRQ(info);
886:   info = xx->GetArray(&tptr2,&nn2);CHKERRQ(info);
887:   info = yy->GetArray(&tptr3,&nn3);CHKERRQ(info);
888:   if (nn1!=nn2 || nn2!=nn3) {TaoFunctionReturn(1);}
889:   for (i=0;i<nn1;i++){ tptr1[i] = a * tptr2[i] + b * tptr3[i]; }
890:   info = this->RestoreArray(&tptr1,&nn1);CHKERRQ(info);
891:   info = xx->RestoreArray(&tptr2,&nn2);CHKERRQ(info);
892:   info = yy->RestoreArray(&tptr3,&nn3);CHKERRQ(info);
893:   TaoFunctionReturn(0);
894: }

898: /*@C
899:    AbsoluteValue - Sets each element of this vector equal to its magnitude.

901:    Input Parameters: none

903:    Level: intermediate
904: @*/
905: int TaoVec::AbsoluteValue(){
906:   int info;
907:   TaoInt nn,i;
908:   TaoScalar *v;

910:   TaoFunctionBegin;
911:   info = this->GetArray(&v,&nn);CHKERRQ(info);
912:   for (i=0;i<nn;i++){ 
913:     v[i]= TaoAbsScalar(v[i]);
914:   }
915:   info = this->RestoreArray(&v,&nn);CHKERRQ(info);
916:   TaoFunctionReturn(0);
917: }

921: /*@C
922:    MinElement - Finds the smallest element of this vector.

924:    Output Parameter:
925: .  val - the smallest value in the vector

927:    Level: intermediate
928: @*/
929: int TaoVec::MinElement(double *val){
930:   int info;
931:   TaoInt nn,i;
932:   TaoScalar dd=TAO_INFINITY,*v;

934:   TaoFunctionBegin;
935:   info = this->GetArray(&v,&nn);CHKERRQ(info);
936:   for (i=0;i<nn;i++){ 
937:     if (v[i]<dd) dd=v[i];
938:   }
939:   info = this->RestoreArray(&v,&nn);CHKERRQ(info);
940:   *val = dd;
941:   TaoFunctionReturn(0);
942: }

946: /*@C
947:    Sign - Replace each element of this vector with -1, 0, or 1, depending on its sign.

949:    Level: intermediate
950: @*/
951: int TaoVec::Sign(){
952:   int info;
953:   TaoInt nn,i;
954:   TaoScalar *v;

956:   TaoFunctionBegin;
957:   info = this->GetArray(&v,&nn);CHKERRQ(info);
958:   for (i=0;i<nn;i++){
959:     if (v[i]<0){ 
960:       v[i]=-1;
961:     } else if (v[i]>0){
962:       v[i]=1.0;
963:     } else {
964:       v[i]=0.0;
965:     }
966:   }
967:   info = this->RestoreArray(&v,&nn);CHKERRQ(info);
968:   TaoFunctionReturn(0);
969: }

973: /*@C
974:    SetReducedVec - Sets the reduced space of the vector that this
975:    vector should represent.  The index set describes which
976:    elements of the vector should be used.  This routine also
977:    copies  the appropriate elements from the full space vector
978:    into the reduced space vector.

980:    Input Parameters:
981: +  vv -  a vector
982: -  ss -  an index set

984:    Level: advanced

986: .seealso TaoVec::CreateReducedVec(),  TaoVec::ReducedCopyFromFull(), TaoVec::ReducedXPY(), TaoSelectSubset()
987: @*/
988: int TaoVec::SetReducedVec(TaoVec* vv,TaoIndexSet* ss){
989:   TaoFunctionBegin;
990:   SETERRQ(56,"Operation not defined");
991:   /* TaoFunctionReturn(1); */
992: }

996: /*@C
997:   ReducedCopyFromFull - Copies the appropriate elements of the vector into
998:   this reduced vector.
999:   
1000:   Input Parameters:
1001: +  ss -  an index set
1002: -  vv -  a full vector

1004:    Level: advanced

1006: .seealso TaoVec::CreateReducedVec(), TaoVec::ReducedXPY()
1007: @*/
1008: int TaoVec::ReducedCopyFromFull(TaoVec* vv, TaoIndexSet* ss){
1009:   TaoFunctionBegin;
1010:   SETERRQ(56,"Operation not defined");
1011:   /* TaoFunctionReturn(1); */
1012: }

1016: /*@C
1017:   ReducedXPY - Adds a reduced vector to the appropriate elements of this vector.
1018:   
1019:   Input Parameters:
1020: +  vv -  the reduced vector
1021: -  ss -  the index set identifying which elements of this vector should be supplemented

1023:    Level: advanced

1025: .seealso TaoVec::CreateReducedVec(), TaoVec::ReducedCopyFromFull()
1026: @*/
1027: int TaoVec::ReducedXPY(TaoVec* vv, TaoIndexSet* ss){
1028:   TaoFunctionBegin;
1029:   SETERRQ(56,"Operation not defined");
1030:   /* TaoFunctionReturn(1); */
1031: }

1035: /*@C
1036:   StepMax - Calculates the largest multiple of a vector that can be added to
1037:   this vector while keeping each element of this vector nonnegative.
1038:   
1039:   Input Parameters:
1040: . vv - the step direction

1042:   Input Parameters:
1043: . step1 - the maximum stepsize

1045:   Note:
1046:   This vector should contain all positive elements.

1048:   Note:
1049:   If there is no maximum steplength, the output scalar may be set
1050:   to TAO_INFINITY.

1052:   Level: advanced
1053: @*/
1054: int TaoVec::StepMax(TaoVec* vv,double *step1){
1055:   int info;
1056:   TaoInt nn1,nn2,i;
1057:   TaoScalar *xx,*dx;
1058:   double stepmax1=TAO_INFINITY;

1060:   TaoFunctionBegin;
1061:   info = this->GetArray(&xx,&nn1);CHKERRQ(info);
1062:   if (this!=vv){
1063:     info = vv->GetArray(&dx,&nn2);CHKERRQ(info);
1064:     if (nn1!=nn2) {TaoFunctionReturn(1);}
1065:     for (i=0;i<nn1;i++){
1066:       if (xx[i] < 0){
1067:         TaoFunctionReturn(1);
1068:       } else if (dx[i]<0){ 
1069:         stepmax1=TaoMin(stepmax1,-xx[i]/dx[i]);
1070:       }
1071:     }
1072:     info = vv->RestoreArray(&dx,&nn2);CHKERRQ(info);
1073:     *step1=stepmax1;
1074:   } else {
1075:     *step1=1.0;
1076:   }
1077:   info = this->RestoreArray(&xx,&nn1);CHKERRQ(info);
1078:   TaoFunctionReturn(0);
1079:   }

1081: /*@C
1082:   StepBoundInfo - Calculates the largest multiple of a vector that can be added to
1083:   this vector while keeping each element of this vector nonnegative.
1084:   
1085:   Input Parameters:
1086: + txl - the lower bounds on this vector
1087: . txu - the upper bounds on this vector
1088: - tdx - the step direction for this vector

1090:   Output Parameters:
1091: + boundmin - the step to closest bound i.e min(a1, ..., an);
1092: . wolfemin - the step to closest bound not equal i.e min(b1, ..., bn);
1093: - boundmax - the step to farthest bound   i.e. max(c1, ..., cn);

1095:   Where:
1096:   if tdx[i] > 0;  ai = (txu[i] - this[i])/tdx[i] ; bi=ai, ci=ai;
1097:   if tdx[i] < 0;  ai = (txl[i] - this[i])/tdx[i] ; bi=ai, ci=ai
1098:   if tdx[i] == 0 && txl[i] < x[i] < txu[i] ; ai=TAO_INFINITY, bi=ai, ci=ai;
1099:   if tdx[i] == 0 && (txl[i] == x[i] || txu[i] == x[i]) ; ai= 0, bi=TAO_INFINITY, ci=0;
1100:  
1101:   Note:
1102:   If there is no maximum steplength, the output scalar may be set
1103:   to TAO_INFINITY.

1105:   Level: advanced
1106: @*/
1107: int TaoVec::StepBoundInfo(TaoVec* XL ,TaoVec* XU, TaoVec*S, double *bmin1,double *bmin2, double *bmax){
1108:   TaoFunctionBegin;
1109:   SETERRQ(56,"Operation not defined");
1110:   /* TaoFunctionReturn(1); */
1111: }

1115: /*@C
1116:   View - Views the contents of the vector.
1117:   
1118:   Input Parameters: none
1119:   
1120:   Level: intermediate
1121: @*/
1122: int TaoVec::View(){
1123:   TaoFunctionBegin;
1124:   SETERRQ(56,"Operation not defined");
1125:   /* TaoFunctionReturn(1); */
1126: }

1130: /*@C
1131:   BoundGradientProjection - Projects this vector according to this definition.
1132:   If XX[i]==XXL[i], then this[i] = min(G[i],0);   
1133:   If XX[i]==XXU[i], then this[i] = max(G[i],0);
1134:   else this[i] = G[i];

1136:   Input Parameters:
1137: . GG,XX,XXL,XXU - the vectors.

1139:   Level: advanced
1140: @*/
1141: int TaoVec::BoundGradientProjection(TaoVec* gg,TaoVec* xxll,TaoVec* xx, TaoVec* xxuu){
1142:   int info;
1143:   TaoInt nn1,nn2,nn3,nn4,nn5,i;
1144:   TaoScalar *xptr,*xlptr,*xuptr,*gptr,*gpptr;

1146:   TaoFunctionBegin;
1147:   info = this->GetArray(&gpptr,&nn1);CHKERRQ(info);
1148:   if (this != gg){
1149:     info = gg->GetArray(&gptr,&nn2);CHKERRQ(info);
1150:   } else {
1151:     gptr=gpptr; nn2=nn1;
1152:   }
1153:   info = xxll->GetArray(&xlptr,&nn3);CHKERRQ(info);
1154:   info = xx->GetArray(&xptr,&nn4);CHKERRQ(info);
1155:   info = xxuu->GetArray(&xuptr,&nn5);CHKERRQ(info);
1156:   if (nn1!=nn2 || nn2!=nn3 || nn3!=nn4 || nn4!=nn5) {TaoFunctionReturn(1);}

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

1160:     gpptr[i] = gptr[i];
1161:     if (gpptr[i]>0 && xptr[i]<=xlptr[i]){
1162:       gpptr[i] = 0;
1163:     } else if (gpptr[i]<0 && xptr[i]>=xuptr[i]){
1164:       gpptr[i] = 0;
1165:     }
1166:   }
1167:   info = this->RestoreArray(&gpptr,&nn1);CHKERRQ(info);
1168:   if (this!=gg){
1169:     info = gg->RestoreArray(&gptr,&nn2);CHKERRQ(info);
1170:   }
1171:   info = xxll->RestoreArray(&xlptr,&nn3);CHKERRQ(info);
1172:   info = xx->RestoreArray(&xptr,&nn4);CHKERRQ(info);
1173:   info = xxuu->RestoreArray(&xuptr,&nn5);CHKERRQ(info);

1175:   TaoFunctionReturn(0);
1176: }

1180: /*@C
1181:    GetArray - Sets a pointer to the first element in the vector array.

1183:    Output Parameters:
1184: +  dptr -  pointer an the array of numbers
1185: -  n -  the length of the array

1187:    Note:
1188:    This operation may not be defined the same for all vector types.

1190:    Level: intermediate

1192: .seealso TaoVec::RestoreArray()
1193: @*/
1194: int TaoVec::GetArray(TaoScalar **dptr, TaoInt *n){
1195:   TaoFunctionBegin;
1196:   SETERRQ(56,"Operation not defined");
1197:   /* TaoFunctionReturn(1); */
1198: }

1202: /*@C
1203:    RestoreArray - Returns a pointer to the first element in the vector array.

1205:    Input Parameters:
1206: +  dptr -  pointer an the array of numbers obtained from TaoVec::GetArray
1207: -  n -  the length of the array

1209:    Note:
1210:    This routine is not used within TAO solvers.  Rather, it is used to
1211:    implement other routines.

1213:    Level: intermediate

1215: .seealso TaoVec::GetArray()
1216: @*/
1217: int TaoVec::RestoreArray(TaoScalar **dptr, TaoInt *n){
1218:   TaoFunctionBegin;
1219:   *dptr=NULL;
1220:   *n=0;
1221:   TaoFunctionReturn(0);
1222: }

1226: /*@C
1227:    CreateIndexSet - Creates an index set that may be used to describe sets of
1228:    elements of this vector.

1230:    Output Parameters:
1231: .  ss -  a new index set

1233:    Level: advanced
1234: @*/
1235: int TaoVec::CreateIndexSet(TaoIndexSet **ss){
1236:   TaoFunctionBegin;
1237:   SETERRQ(56,"Operation not defined");
1238:   /* TaoFunctionReturn(1); */
1239: }

1241: inline static double fischer(double a, double b) 
1242: {
1243:    // Method suggested by Bob Vanderbei
1244:    if (a + b <= 0) {
1245:      return sqrt(a*a + b*b) - (a + b);
1246:    }
1247:    return -2.0*a*b / (sqrt(a*a + b*b) + (a + b));
1248: }

1252: /*@C
1253:    Fischer - Evaluates the Fischer-Burmeister function for complementarity 
1254:    problems.
1255:    
1256:    Input Parameters:
1257: +  xx - current point
1258: .  ff - function evaluated at x
1259: .  ll - lower bounds 
1260: -  uu - upper bounds

1262:    Notes: 
1263:    The Fischer-Burmeister function is defined as
1264: $        phi(a,b) := sqrt(a*a + b*b) - a - b
1265:    and is used reformulate a complementarity problem as a semismooth
1266:    system of equations.

1268:    The result of this function is done by cases:
1269: +  l[i] == -infinity, u[i] == infinity  -- fb[i] = -f[i]
1270: .  l[i] == -infinity, u[i] finite       -- fb[i] = phi(u[i]-x[i], -f[i])
1271: .  l[i] finite,       u[i] == infinity  -- fb[i] = phi(x[i]-l[i],  f[i])
1272: .  l[i] finite < u[i] finite -- fb[i] = phi(x[i]-l[i], phi(u[i]-x[i], -f[u]))
1273: -  otherwise l[i] == u[i] -- fb[i] = l[i] - x[i]

1275:    Level: advanced

1277: .seealso TaoMat::D_Fischer()
1278: @*/
1279: int TaoVec::Fischer(TaoVec* xx, TaoVec* ff, TaoVec* ll, TaoVec* uu)
1280: {
1281:   int info;
1282:   TaoInt nn1,nn2,nn3,nn4,nn5,i;
1283:   TaoScalar *v,*x,*f,*l,*u;

1285:   TaoFunctionBegin;
1286:   info = this->GetArray(&v,&nn1); CHKERRQ(info);
1287:   info = xx->GetArray(&x,&nn2); CHKERRQ(info);
1288:   info = ff->GetArray(&f,&nn3); CHKERRQ(info);
1289:   info = ll->GetArray(&l,&nn4); CHKERRQ(info);
1290:   info = uu->GetArray(&u,&nn5); CHKERRQ(info);
1291:   if (nn1!=nn2 || nn2!=nn3 || nn3!=nn4 || nn4!=nn5) { 
1292:     TaoFunctionReturn(1); 
1293:   }

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

1297:     if ((l[i] <= -TAO_INFINITY) && (u[i] >= TAO_INFINITY)) {
1298:       v[i] = -f[i];
1299:     } 
1300:     else if (l[i] <= -TAO_INFINITY) {
1301:       v[i] = -fischer(u[i] - x[i], -f[i]);
1302:     } 
1303:     else if (u[i] >=  TAO_INFINITY) {
1304:       v[i] =  fischer(x[i] - l[i],  f[i]);
1305:     } 
1306:     else if (l[i] == u[i]) {
1307:       v[i] = l[i] - x[i];
1308:     } 
1309:     else {
1310:       v[i] =  fischer(u[i] - x[i], -f[i]);
1311:       v[i] =  fischer(x[i] - l[i],  v[i]);
1312:     }

1314:   }

1316:   info = this->RestoreArray(&v,&nn1);CHKERRQ(info);
1317:   info = xx->RestoreArray(&x,&nn2);CHKERRQ(info);
1318:   info = ff->RestoreArray(&f,&nn3);CHKERRQ(info);
1319:   info = ll->RestoreArray(&l,&nn4);CHKERRQ(info);
1320:   info = uu->RestoreArray(&u,&nn5);CHKERRQ(info);

1322:   TaoFunctionReturn(0);
1323: }

1325: inline static double sfischer(double a, double b, double c)
1326: {
1327:    // Method suggested by Bob Vanderbei
1328:    if (a + b <= 0) {
1329:      return sqrt(a*a + b*b + 2.0*c*c) - (a + b);
1330:    }
1331:    return 2.0*(c*c - a*b) / (sqrt(a*a + b*b + 2.0*c*c) + (a + b));
1332: }

1336: /*@C
1337:    SFischer - Evaluates the Smoothed Fischer-Burmeister function for
1338:    complementarity problems.

1340:    Input Parameters:
1341: +  xx - current point
1342: .  ff - function evaluated at x
1343: .  ll - lower bounds
1344: .  uu - upper bounds
1345: -  mu - smoothing parameter

1347:    Notes:
1348:    The Smoothed Fischer-Burmeister function is defined as
1349: $        phi(a,b) := sqrt(a*a + b*b + 2*mu*mu) - a - b
1350:    and is used reformulate a complementarity problem as a semismooth
1351:    system of equations.

1353:    The result of this function is done by cases:
1354: +  l[i] == -infinity, u[i] == infinity  -- fb[i] = -f[i] - 2*mu*x[i]
1355: .  l[i] == -infinity, u[i] finite       -- fb[i] = phi(u[i]-x[i], -f[i], mu)
1356: .  l[i] finite,       u[i] == infinity  -- fb[i] = phi(x[i]-l[i],  f[i], mu)
1357: .  l[i] finite < u[i] finite -- fb[i] = phi(x[i]-l[i], phi(u[i]-x[i], -f[u], mu), mu)
1358: -  otherwise l[i] == u[i] -- fb[i] = l[i] - x[i]

1360:    Level: advanced

1362: .seealso TaoMat::SD_Fischer()
1363: @*/
1364: int TaoVec::SFischer(TaoVec* xx, TaoVec* ff, TaoVec* ll, TaoVec* uu, double mu)
1365: {

1367:   int info;
1368:   TaoInt nn1, nn2, nn3, nn4, nn5,i;
1369:   TaoScalar *v, *x, *f, *l, *u;

1371:   TaoFunctionBegin;

1373:   if ((mu >= -TAO_EPSILON) && (mu <= TAO_EPSILON)) {
1374:     Fischer(xx, ff, ll, uu);
1375:   }
1376:   else {
1377:     info = this->GetArray(&v, &nn1); CHKERRQ(info);
1378:     info = xx->GetArray(&x, &nn2); CHKERRQ(info);
1379:     info = ff->GetArray(&f, &nn3); CHKERRQ(info);
1380:     info = ll->GetArray(&l, &nn4); CHKERRQ(info);
1381:     info = uu->GetArray(&u, &nn5); CHKERRQ(info);

1383:     if (nn1!=nn2 || nn2!=nn3 || nn3!=nn4 || nn4!=nn5) {
1384:       TaoFunctionReturn(1);
1385:     }

1387:     for (i = 0; i < nn1; ++i) {
1388:       if ((l[i] <= -TAO_INFINITY) && (u[i] >= TAO_INFINITY)) {
1389:         v[i] = -f[i] - mu*x[i];
1390:       } 
1391:       else if (l[i] <= -TAO_INFINITY) {
1392:         v[i] = -sfischer(u[i] - x[i], -f[i], mu);
1393:       } 
1394:       else if (u[i] >=  TAO_INFINITY) {
1395:         v[i] =  sfischer(x[i] - l[i],  f[i], mu);
1396:       } 
1397:       else if (l[i] == u[i]) {
1398:         v[i] = l[i] - x[i];
1399:       } 
1400:       else {
1401:         v[i] =  sfischer(u[i] - x[i], -f[i], mu);
1402:         v[i] =  sfischer(x[i] - l[i],  v[i], mu);
1403:       }
1404:     }

1406:     info = this->RestoreArray(&v, &nn1); CHKERRQ(info);
1407:     info = xx->RestoreArray(&x, &nn2); CHKERRQ(info);
1408:     info = ff->RestoreArray(&f, &nn3); CHKERRQ(info);
1409:     info = ll->RestoreArray(&l, &nn4); CHKERRQ(info);
1410:     info = uu->RestoreArray(&u, &nn5); CHKERRQ(info);
1411:   }
1412:   TaoFunctionReturn(0);
1413: }