ScalarVectors.c 3.67 KB
Newer Older
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
#include <stdio.h>
#include <stdlib.h>
#include <math.h>
#include <string.h>
#include <ScalarVectors.h>

/*********************************************************************************/

void CreateInts (int **vint, int dim) {
	if ((*vint = (int *) malloc (sizeof(int)*dim)) == NULL)
		{ printf ("Memory Error (CreateInts(%d))\n", dim); exit (1); }
}

void RemoveInts (int **vint) { 
	if (*vint != NULL) { free (*vint); *vint = NULL; }  
}

void InitInts (int *vint, int dim, int frst, int incr) {
	int i, *p1 = vint, num = frst;

	for (i=0; i<dim; i++) 
		{ *(p1++) = num; num += incr; }
}
		
void CopyInts (int *src, int *dst, int dim) { 
	memmove (dst, src, sizeof(int) * dim);
}

void CopyShiftInts (int *src, int *dst, int dim, int shft) {
	int i, *p1 = src, *p2 = dst;

	if (shft == 0)
		CopyInts (src, dst, dim);
	else
		for (i=0; i<dim; i++)
			*(p2++) = *(p1++) + shft;
}
	
void TransformLengthtoHeader (int *vint, int dim) {
	int i, *pi = vint; 

	for (i=0; i<dim; i++) { *(pi+1) += *pi; pi++; }
}

void TransformHeadertoLength (int *vint, int dim) {
	int i, *pi = vint+dim; 

	for (i=dim; i>0; i--) { *(pi) -= *(pi-1); pi--; }
}

void ComputeHeaderfromLength (int *len, int *head, int dim) {
	int i, *pi1 = len, *pi2 = head; 

	for (i=0; i<dim; i++) { *(pi2+1) = (*pi2) +(*(pi1++)); pi2++; }
}

void ComputeLengthfromHeader (int *head, int *len, int dim) {
	int i, *pi1 = head, *pi2 = len; 

	for (i=0; i<dim; i++) { *(pi2++) = (*(pi1+1)) -(*pi1); pi1++; }
}

int AddInts (int *vint, int dim) {
	int i, *pi = vint, aux = 0;

	for (i=0; i<dim; i++) { 
		aux += *pi; pi++; 
	}

	return aux;
}

// The permutation defined by perm is applied on vec, whose size is dim. 
void PermuteInts (int *vec, int *perm, int dim) {
  int i, *pi = vec;

  for (i=0; i<dim; i++) { *pi = perm[*pi]; pi++; }
}

// Apply the inverse of perm, and store it on iperm, whose size is dim. 
void ComputeInvPermutation (int *perm, int *iperm, int dim) {
  int i, *pi1 = perm;

  for (i=0; i<dim; i++) { iperm[*(pi1++)] = i; }
}

// Scale by scal the elements of vint, whose size is dim. 
void ScaleInts (int *vint, int scal, int dim) {
  int i;
  int *pi = vint;

  for (i=0; i<dim; i++)
    *(pi++) *= scal;
}

/*********************************************************************************/

void CreateDoubles (double **vdbl, int dim) {
	if ((*vdbl = (double *) malloc (sizeof(double)*dim)) == NULL)
		{ printf ("Memory Error (CreateDoubles(%d))\n", dim); exit (1); }
}

void RemoveDoubles (double **vdbl) { 
	if (*vdbl != NULL) { free (*vdbl); *vdbl = NULL; }
}

void InitDoubles (double *vdbl, int dim, double frst, double incr) {
	int i; 
	double *pd = vdbl, num = frst;

	for (i=0; i<dim; i++) 
		{ *(pd++) = num; num += incr; }
}
		
void InitRandDoubles (double *vdbl, int dim, double frst, double last) {
	int i; 
	double *pd = vdbl, size = last - frst;

	for (i=0; i<dim; i++) 
		{ *(pd++) = frst + (size * (rand() / (RAND_MAX + 1.0))); }
}
		
void CopyDoubles (double *src, double *dst, int dim) { 
	memmove (dst, src, sizeof(double) * dim);
}

void ScaleDoubles (double *vdbl, double scal, int dim) {
	int i; 
	double *pd = vdbl;

	for (i=0; i<dim; i++) 
		*(pd++) *= scal;
}

double DotDoubles (double *vdbl1, double *vdbl2, int dim) {
	int i; 
	double *pd1 = vdbl1, *pd2 = vdbl2, res = 0.0;

	for (i=0; i<dim; i++) 
		res += (*(pd2++)) * (*(pd1++));

	return res;
}

void VvecDoubles (double alfa, double *src1, double *src2, double beta, double *dst, int dim) {
    int i;

    for (i = 0; i < dim; i++) {
        //dst[i] = (beta * dst[i]) + (alfa * src1[i] * src2[i]); 
        double tmp = alfa * src1[i] * src2[i];
        dst[i] = fma(beta, dst[i], tmp);
    }
}


/*********************************************************************************/