# Basic linear computations

Last Update: Feb. 9, 2006

## MPFR

### Sample code

basic_linear_mpf.c

```#include <stdio.h>
#include <math.h>

/* Header file of BNCpack */
#include "bnc.h"

/* Basic linear computations */

#define DIM 5 /* Dimension of Vectors and a matrix */

/* Main function */
int main()
{
long int i, j;
mpf_t ret, tmp;

/* multiple-precision real vectors */
MPFVector mpfva, mpfvb, mpfvc;
/* multiple-precision real square matrix */
MPFMatrix mpfma;

/* Set default prec in decimal */
set_bnc_default_prec_decimal(50);

/* Initialize */
mpf_init(ret);
mpf_init(tmp);

/* Initialize vectors */
mpfva = init_mpfvector(DIM);
mpfvb = init_mpfvector(DIM);
mpfvc = init_mpfvector(DIM);

/* Initialize a matrix */
mpfma = init_mpfmatrix(DIM, DIM);

/* Set elements of vectors */
for(i = 0; i < DIM; i++)
{
set_mpfvector_i_ui(mpfva, i, i + 1);
set_mpfvector_i_ui(mpfvb, i, DIM - i);
}

/* Set elements of a matrix */
for(i = 0; i < DIM; i++)
{
for(j = 0; j < DIM; j++)
set_mpfmatrix_ij_ui(mpfma, i, j, i + j + 1);
}

/* print a matrix */
printf("mpfma = \n"); print_mpfmatrix(mpfma); printf("\n");

/* print vectors */
printf("mpfva = \n"); print_mpfvector(mpfva); printf("\n");
printf("mpfvb = \n"); print_mpfvector(mpfvb); printf("\n");

printf("mpfva + mpfvb = \n"); print_mpfvector(mpfvc); printf("\n");

/* subtraction of vectors */
sub_mpfvector(mpfvc, mpfva, mpfvb);
printf("mpfva - mpfvb = \n"); print_mpfvector(mpfvc); printf("\n");

/* a scalar times a vector */
mpf_sqrt_ui(tmp, 2UL);
cmul_mpfvector(mpfvc, tmp, mpfva);
printf("sqrt(2) * mpfva = \n"); print_mpfvector(mpfvc); printf("\n");

/* inner product of vectors */
ip_mpfvector(ret, mpfva, mpfvb);
printf("(mpfva, mpfvb)= "); mpf_out_str(stdout, 10, 0, ret); printf("\n");

/* 2-norm of a vector */
norm2_mpfvector(ret, mpfva);
printf("||mpfva||_2   = "); mpf_out_str(stdout, 10, 0, ret); printf("\n");

/* 1-norm of a vector */
norm1_mpfvector(ret, mpfva);
printf("||mpfva||_1   = "); mpf_out_str(stdout, 10, 0, ret); printf("\n");

/* inf-norm of a vector */
normi_mpfvector(ret, mpfva);
printf("||mpfva||_inf = "); mpf_out_str(stdout, 10, 0, ret); printf("\n\n");

/* matrix-vector product */
mul_mpfmatrix_mpfvec(mpfvc, mpfma, mpfva);
printf("mpfma * mpfva = \n"); print_mpfvector(mpfvc); printf("\n");

/* Clear */
mpf_clear(ret);
mpf_clear(tmp);

/* Clear vectors */
free_mpfvector(mpfva);
free_mpfvector(mpfvb);
free_mpfvector(mpfvc);

/* Clear a matrix */
free_mpfmatrix(mpfma);

return 0;
}```

### Example execution

```\$ ./basic_linear_mpf
-------------------------------------------------------------------------------
BNC Default Precision    : 167 bits(50.3 decimal digits)
BNC Default Rounding Mode: Round to Nearest
-------------------------------------------------------------------------------
mpfma =
0 1.000000000000000000000000000000000000000000000000000 2.000000000000000000000000000000000000000000000000000 3.000000000000000000000000000000000000000000000000000 4.000000000000000000000000000000000000000000000000000 5.000000000000000000000000000000000000000000000000000
1 2.000000000000000000000000000000000000000000000000000 3.000000000000000000000000000000000000000000000000000 4.000000000000000000000000000000000000000000000000000 5.000000000000000000000000000000000000000000000000000 6.000000000000000000000000000000000000000000000000000
2 3.000000000000000000000000000000000000000000000000000 4.000000000000000000000000000000000000000000000000000 5.000000000000000000000000000000000000000000000000000 6.000000000000000000000000000000000000000000000000000 7.000000000000000000000000000000000000000000000000000
3 4.000000000000000000000000000000000000000000000000000 5.000000000000000000000000000000000000000000000000000 6.000000000000000000000000000000000000000000000000000 7.000000000000000000000000000000000000000000000000000 8.000000000000000000000000000000000000000000000000000
4 5.000000000000000000000000000000000000000000000000000 6.000000000000000000000000000000000000000000000000000 7.000000000000000000000000000000000000000000000000000 8.000000000000000000000000000000000000000000000000000 9.000000000000000000000000000000000000000000000000000

mpfva =
0 1.000000000000000000000000000000000000000000000000000
1 2.000000000000000000000000000000000000000000000000000
2 3.000000000000000000000000000000000000000000000000000
3 4.000000000000000000000000000000000000000000000000000
4 5.000000000000000000000000000000000000000000000000000

mpfvb =
0 5.000000000000000000000000000000000000000000000000000
1 4.000000000000000000000000000000000000000000000000000
2 3.000000000000000000000000000000000000000000000000000
3 2.000000000000000000000000000000000000000000000000000
4 1.000000000000000000000000000000000000000000000000000

mpfva + mpfvb =
0 6.000000000000000000000000000000000000000000000000000
1 6.000000000000000000000000000000000000000000000000000
2 6.000000000000000000000000000000000000000000000000000
3 6.000000000000000000000000000000000000000000000000000
4 6.000000000000000000000000000000000000000000000000000

mpfva - mpfvb =
0 -4.000000000000000000000000000000000000000000000000000
1 -2.000000000000000000000000000000000000000000000000000
2 0
3 2.000000000000000000000000000000000000000000000000000
4 4.000000000000000000000000000000000000000000000000000

sqrt(2) * mpfva =
0 1.414213562373095048801688724209698078569671875376952
1 2.828427124746190097603377448419396157139343750753904
2 4.242640687119285146405066172629094235709015626130867
3 5.656854249492380195206754896838792314278687501507809
4 7.071067811865475244008443621048490392848359376884750

(mpfva, mpfvb)= 3.500000000000000000000000000000000000000000000000000e1
||mpfva||_2   = 7.416198487095662948711397440800713060979904319097515
||mpfva||_1   = 1.500000000000000000000000000000000000000000000000000e1
||mpfva||_inf = 5.000000000000000000000000000000000000000000000000000

mpfma * mpfva =
0 5.500000000000000000000000000000000000000000000000000e1
1 7.000000000000000000000000000000000000000000000000000e1
2 8.500000000000000000000000000000000000000000000000000e1
3 1.000000000000000000000000000000000000000000000000000e2
4 1.150000000000000000000000000000000000000000000000000e2```

## IEEE754 double

### Sample code

basic_linear_d.c

```#include <stdio.h>
#include <math.h>

/* Header file of BNCpack */
#include "bnc.h"

/* Basic linear computations */

#define DIM 5 /* Dimension of Vectors and a matrix */

/* Main function */
int main()
{
long int i, j;

/* IEEE754 double real vectors */
DVector dva, dvb, dvc;
/* IEEE754 double real square matrix */
DMatrix dma;

/* Initialize vectors */
dva = init_dvector(DIM);
dvb = init_dvector(DIM);
dvc = init_dvector(DIM);

/* Initialize a matrix */
dma = init_dmatrix(DIM, DIM);

/* Set elements of vectors */
for(i = 0; i < DIM; i++)
{
set_dvector_i(dva, i, (double)(i + 1));
set_dvector_i(dvb, i, (double)(DIM - i));
}

/* Set elements of a matrix */
for(i = 0; i < DIM; i++)
{
for(j = 0; j < DIM; j++)
set_dmatrix_ij(dma, i, j, (double)(i + j + 1));
}

/* print a matrix */
printf("dma = \n"); print_dmatrix(dma); printf("\n");

/* print vectors */
printf("dva = \n"); print_dvector(dva); printf("\n");
printf("dvb = \n"); print_dvector(dvb); printf("\n");

printf("dva + dvb = \n"); print_dvector(dvc); printf("\n");

/* subtraction of vectors */
sub_dvector(dvc, dva, dvb);
printf("dva - dvb = \n"); print_dvector(dvc); printf("\n");

/* a scalar times a vector */
cmul_dvector(dvc, sqrt(2.0), dva);
printf("sqrt(2) * dva = \n"); print_dvector(dvc); printf("\n");

/* inner product of vectors */
printf("(dva, dvb)  = %25.17e\n", ip_dvector(dva, dvb));

/* 2-norm of a vector */
printf("||dva||_2   = %25.17e\n", norm2_dvector(dva));

/* 1-norm of a vector */
printf("||dva||_1   = %25.17e\n", norm1_dvector(dva));

/* inf-norm of a vector */
printf("||dva||_inf = %25.17e\n", normi_dvector(dva)); printf("\n");

/* matrix-vector product */
mul_dmatrix_dvec(dvc, dma, dva);
printf("dma * dva = \n"); print_dvector(dvc); printf("\n");

/* Clear vectors */
free_dvector(dva);
free_dvector(dvb);
free_dvector(dvc);

/* Clear a matrix */
free_dmatrix(dma);

return 0;
}```

### Example execution

```\$ ./basic_linear_d
dma =
0   1.00000000000000000e+00   2.00000000000000000e+00   3.00000000000000000e+00   4.00000000000000000e+00   5.00000000000000000e+00
1   2.00000000000000000e+00   3.00000000000000000e+00   4.00000000000000000e+00   5.00000000000000000e+00   6.00000000000000000e+00
2   3.00000000000000000e+00   4.00000000000000000e+00   5.00000000000000000e+00   6.00000000000000000e+00   7.00000000000000000e+00
3   4.00000000000000000e+00   5.00000000000000000e+00   6.00000000000000000e+00   7.00000000000000000e+00   8.00000000000000000e+00
4   5.00000000000000000e+00   6.00000000000000000e+00   7.00000000000000000e+00   8.00000000000000000e+00   9.00000000000000000e+00

dva =
0   1.00000000000000000e+00
1   2.00000000000000000e+00
2   3.00000000000000000e+00
3   4.00000000000000000e+00
4   5.00000000000000000e+00

dvb =
0   5.00000000000000000e+00
1   4.00000000000000000e+00
2   3.00000000000000000e+00
3   2.00000000000000000e+00
4   1.00000000000000000e+00

dva + dvb =
0   6.00000000000000000e+00
1   6.00000000000000000e+00
2   6.00000000000000000e+00
3   6.00000000000000000e+00
4   6.00000000000000000e+00

dva - dvb =
0  -4.00000000000000000e+00
1  -2.00000000000000000e+00
2   0.00000000000000000e+00
3   2.00000000000000000e+00
4   4.00000000000000000e+00

sqrt(2) * dva =
0   1.41421356237309515e+00
1   2.82842712474619029e+00
2   4.24264068711928566e+00
3   5.65685424949238058e+00
4   7.07106781186547551e+00

(dva, dvb)  =   3.50000000000000000e+01
||dva||_2   =   7.41619848709566298e+00
||dva||_1   =   1.50000000000000000e+01
||dva||_inf =   5.00000000000000000e+00

dma * dva =
0   5.50000000000000000e+01
1   7.00000000000000000e+01
2   8.50000000000000000e+01
3   1.00000000000000000e+02
4   1.15000000000000000e+02```

<- Go back