/[escript]/trunk/escript/src/DataArrayView.h
ViewVC logotype

Diff of /trunk/escript/src/DataArrayView.h

Parent Directory Parent Directory | Revision Log Revision Log | View Patch Patch

revision 576 by gross, Fri Mar 3 08:28:42 2006 UTC revision 580 by gross, Wed Mar 8 05:45:51 2006 UTC
# Line 773  class DataArrayView { Line 773  class DataArrayView {
773       \brief       \brief
774       solves a local eigenvalue problem       solves a local eigenvalue problem
775    
776       \param in - Input - The left hand side.       \param in - Input - matrix
777         \param inOffset - Input - offset into in
778       \param ev - Output - The eigenvalues       \param ev - Output - The eigenvalues
779         \param inOffset - Input - offset into ev
780    */    */
781      static
782    inline    inline
783    void    void
784    DataArrayView::eigenvalues(const DataArrayView& in,    DataArrayView::eigenvalues(DataArrayView& in,
785                               DataArrayView& ev);                               ValueType::size_type inOffset,
786                                 DataArrayView& ev,
787                                 ValueType::size_type evOffset)
788     {
789       double in00,in10,in20,in01,in11,in21,in02,in12,in22;
790       double ev0,ev1,ev2;
791       int s=in.getShape()[0];
792       if (s==1) {
793          in00=(*(in.m_data))[inOffset+in.index(0,0)];
794          eigenvalues1(in00,&ev0);
795          (*(ev.m_data))[evOffset+ev.index(0)]=ev0;
796    
797       } else  if (s==2) {
798          in00=(*(in.m_data))[inOffset+in.index(0,0)];
799          in10=(*(in.m_data))[inOffset+in.index(1,0)];
800          in01=(*(in.m_data))[inOffset+in.index(0,1)];
801          in11=(*(in.m_data))[inOffset+in.index(1,1)];
802          eigenvalues2(in00,(in01+in10)/2.,in11,&ev0,&ev1);
803          (*(ev.m_data))[evOffset+ev.index(0)]=ev0;
804          (*(ev.m_data))[evOffset+ev.index(1)]=ev1;
805    
806       } else  if (s==3) {
807          in00=(*(in.m_data))[inOffset+in.index(0,0)];
808          in10=(*(in.m_data))[inOffset+in.index(1,0)];
809          in20=(*(in.m_data))[inOffset+in.index(2,0)];
810          in01=(*(in.m_data))[inOffset+in.index(0,1)];
811          in11=(*(in.m_data))[inOffset+in.index(1,1)];
812          in21=(*(in.m_data))[inOffset+in.index(2,1)];
813          in02=(*(in.m_data))[inOffset+in.index(0,2)];
814          in12=(*(in.m_data))[inOffset+in.index(1,2)];
815          in22=(*(in.m_data))[inOffset+in.index(2,2)];
816          eigenvalues3(in00,(in01+in10)/2.,(in02+in20)/2.,in11,(in21+in12)/2.,in22,
817                     &ev0,&ev1,&ev2);
818          (*(ev.m_data))[evOffset+ev.index(0)]=ev0;
819          (*(ev.m_data))[evOffset+ev.index(1)]=ev1;
820          (*(ev.m_data))[evOffset+ev.index(2)]=ev2;
821    
822       }
823     }
824    
825    /**    /**
826       \brief       \brief
827         solves a local eigenvalue problem
828    
829         \param in - Input - matrix
830         \param inOffset - Input - offset into in
831         \param ev - Output - The eigenvalues
832         \param evOffset - Input - offset into ev
833         \param V - Output - The eigenvectors
834         \param VOffset - Input - offset into V
835         \param tol - Input - eigenvalues with relative difference tol are treated as equal
836      */
837      static
838      inline
839      void
840      DataArrayView::eigenvalues_and_eigenvectors(DataArrayView& in,
841                                                  ValueType::size_type inOffset,
842                                                  DataArrayView& ev,
843                                                  ValueType::size_type evOffset,
844                                                  DataArrayView& V,
845                                                  ValueType::size_type VOffset,
846                                                  const double tol=1.e-13)
847      {
848       double ev0,ev1,ev2;
849       int s=in.getShape()[0];
850       if (s==1) {
851          eigenvalues1(in(0,0),&ev0);
852          ev(0)=ev0;
853    
854       } else  if (s==2) {
855          eigenvalues2(in(0,0),(in(0,1)+in(1,0))/2.,in(1,1),
856                       &ev0,&ev1);
857          ev(0)=ev0;
858          ev(1)=ev1;
859    
860       } else  if (s==3) {
861          eigenvalues3(in(0,0),(in(0,1)+in(1,0))/2.,(in(0,2)+in(2,0))/2.,in(1,1),(in(2,1)+in(1,2))/2.,in(2,2),
862                     &ev0,&ev1,&ev2);
863          ev(0)=ev0;
864          ev(1)=ev1;
865          ev(2)=ev2;
866    
867       }
868     }
869     /**
870         \brief
871       Perform a matrix multiply of the given views.       Perform a matrix multiply of the given views.
872       This is purely a utility method and has no bearing on this view.       This is purely a utility method and has no bearing on this view.
873    

Legend:
Removed from v.576  
changed lines
  Added in v.580

  ViewVC Help
Powered by ViewVC 1.1.26