matrix_vector.cpp

Go to the documentation of this file.

#include "matrix_vector.h"

#include <math.h>


namespace mumech {

  
double vektor :: norma () const
{
  long i;
  double norma=0.0;
  
  for(i=0;i<p;i++)
    norma+=v[i]*v[i];
  
  return sqrt(norma);
}

void matice :: inverze_3x3 (matice &d)
{
  if(r==3 && s==3){
    int i,j;
    double Det=d.Determinant_3x3();
    
    for(i=0;i<3;i++)
      for(j=0;j<3;j++){
    g(i,j)=1.0*pow(-1.0,i+j)*d.Determinant_3x3_doplnek(j,i)/Det;
      }
  }
  else{
    printf("\nCHYBA - spatny rozmer matice pro pouziti funkce inverze_3x3\n");
  }
}

void enter()
{
  fflush(stdin);scanf("%*c");
}

void e (int i)
{
  printf("\n%d",i);fflush(stdin);scanf("%*c");printf("\n");
}


int MaticexVektor(vektor &res, matice &mat,vektor &vek)
{
  if (mat.s != vek.p) return 0;
  int u,r;
  res.vynulovat();
  for(u=0;u<mat.r;u++)for(r=0;r<mat.s;r++)res.v[u]+=mat.g(u,r)*vek.v[r];
  return 1; 
}

int RozdilVektoru(vektor &a,vektor &b,vektor &c)
{
  int i;
  if(a.p==b.p && b.p==c.p){
    for(i=0;i<a.p;i++)c.v[i]=a.v[i]-b.v[i];
    return 1;
  }
  return 0;
}

void PrintVektor(vektor &vekt)
{
  int u;
  for(u=0;u<vekt.p;u++){
    printf("%d) %4.4f\n",u,vekt.v[u]);
  } 
  printf("\n"); 
}

void PrintMatrix(matice &mat,FILE *&soubor)
{
  int u,v;
  for(u=0;u<mat.r;u++){
    for(v=0;v<mat.s;v++){fprintf(soubor,"%5.4f ",mat.g(u,v));}
    printf("\n");
  } 
  printf("\n"); 
}

void rotate_vector (vektor &V, double angle)
{
  if (angle!=0.0) {
    vektor zaloha(3);
                matice A(3,3);
    double c,s;
    
    c=cos(angle);
    s=sin(angle);
    A.g(0,0)=A.g(1,1)=c*c;
    A.g(1,0)=A.g(0,1)=s*s;
    A.g(2,2)=c*c-s*s;
    A.g(0,2)=c*s;
    A.g(1,2)=-c*s;
    A.g(2,0)=-2.0*c*s;
    A.g(2,1)=2.0*c*s;
    V.g(2)*=2.0;
    zaloha.rovna_se(V);
    MaticexVektor(V,A,zaloha);
    V.g(2)/=2.0;
  }
}
void change_radky(vektor &a, vektor &b) {
                vektor zaloha(a.p);
                zaloha.rovna_se(a);
                a.rovna_se(b);
                b.rovna_se(zaloha);
}
void gaussovaEliminace(matice &m, vektor &right) {
                if (m.r != m.s || right.p != m.r) { printf("Matice a vektor pro Gaussovu eliminaci nemaji odpovidajici rozmery."); }
                int i, j, k, radku = m.r;
                long double nasobitel, delitel;
                bool found;
                vektor *radek;
                vektor zaloha(radku);
                radek = new vektor[radku];
                for (i = 0; i < radku; i++)radek[i].define(radku + 1);
                for (i = 0; i < radku; i++)for (j = 0; j < radku; j++)radek[i].g(j) = 1.0*m.g(i, j);
                for (i = 0; i < radku; i++)radek[i].g(radku) = 1.0*right.g(i);
                for (i = 0; i < radku - 1; i++) {
                                if (radek[i].g(i) == 0.0) {
                                                found = false;
                                                k = i + 1;
                                                while (found == false && k < radku) {
                                                                if (radek[k].g(i) != 0.0) {
                                                                                change_radky(radek[i], radek[k]);
                                                                                found = true;
                                                                }
                                                                k++;
                                                }
                                                if (k == radku && found == false) { printf("Error pri hledani nenuloveho prvku na diagonale."); enter();}
                                }
                                for (j = i + 1; j < radku; j++) {
                                                if (radek[j].g(i) != 0.0) {
                                                                nasobitel = 1.0*radek[i].g(i);
                                                                delitel = 1.0*radek[j].g(i);
                                                                radek[j].krat(1.0*nasobitel);
                                                                radek[j].deleno(1.0*delitel);
                                                                radek[j].odecist(radek[i]);
                                                }
                                }
                }
                //pod diagonalou nuly
                for (i = 0; i < radku; i++) {
                                delitel = 1.0*radek[i].g(i);
                                radek[i].deleno(1.0*delitel);
                }
                //na diagonale jednicky
                for (i = radku - 2; i > -1; i--) {
                                for (j = i + 1; j < radku; j++)radek[i].g(radku) -= radek[j].g(radku)*radek[i].g(j);
                                for (j = i + 1; j < radku; j++)radek[i].g(j) = 0.0;
                }

                //pouze jednicky na diagonale

                /*for (i = 0; i < radku; i++){
                for (j = 0; j < radku + 1; j++)printf("%5.3lf ",radek[i].g(j));
                printf("\n");
                }
                enter();*/


                for (i = 0; i < radku; i++)right.g(i) = radek[i].g(radku);
};

} // namespace mumech