Skip to content

Commit

Permalink
Merge branch 'use_vector_for_tvector'
Browse files Browse the repository at this point in the history
  • Loading branch information
bastien-roucaries committed Jan 23, 2015
2 parents 13ae936 + 28b97ab commit d983fc5
Show file tree
Hide file tree
Showing 11 changed files with 219 additions and 269 deletions.
75 changes: 43 additions & 32 deletions qucs-core/src/hbsolver.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -26,6 +26,8 @@
# include <config.h>
#endif

#include<algorithm>

#include <stdio.h>

#include "object.h"
Expand Down Expand Up @@ -330,34 +332,34 @@ bool hbsolver::isExcitation (circuit * c) {

// Expands the frequency array using the given frequency and the order.
void hbsolver::expandFrequencies (nr_double_t f, int n) {
tvector<nr_double_t> nfreqs = negfreqs;
tvector<nr_double_t> pfreqs = posfreqs;
int i, k, len = nfreqs.getSize ();
auto nfreqs = negfreqs;
auto pfreqs = posfreqs;
int i, k, len = nfreqs.size ();
negfreqs.clear ();
posfreqs.clear ();
if (len > 0) {
// frequency expansion for full frequency sets
for (i = 0; i <= n + 1; i++) {
for (k = 0; k < len; k++) {
negfreqs.add (i * f + nfreqs.get (k));
negfreqs.push_back (i * f + nfreqs[k]);
}
}
for (i = -n; i < 0; i++) {
for (k = 0; k < len; k++) {
negfreqs.add (i * f + nfreqs.get (k));
negfreqs.push_back (i * f + nfreqs[k]);
}
}
for (i = 0; i <= 2 * n + 1; i++) {
for (k = 0; k < len; k++) {
posfreqs.add (i * f + pfreqs.get (k));
posfreqs.push_back (i * f + pfreqs[k]);
}
}
}
else {
// first frequency
for (i = 0; i <= n + 1; i++) negfreqs.add (i * f);
for (i = -n; i < 0; i++) negfreqs.add (i * f);
for (i = 0; i <= 2 * n + 1; i++) posfreqs.add (i * f);
for (i = 0; i <= n + 1; i++) negfreqs.push_back (i * f);
for (i = -n; i < 0; i++) negfreqs.push_back (i * f);
for (i = 0; i <= 2 * n + 1; i++) posfreqs.push_back (i * f);
}
}

Expand Down Expand Up @@ -388,25 +390,32 @@ void hbsolver::collectFrequencies (void) {
for (auto * c : excitations) {
if (c->getType () != CIR_VDC) { // no extra DC sources
if ((f = c->getPropertyDouble ("f")) != 0.0) {
if (!dfreqs.contains (f)) { // no double frequencies
dfreqs.add (f);
const auto epsilon = std::numeric_limits<nr_double_t>::epsilon();
auto found =
std::find_if(dfreqs.cbegin(),dfreqs.cend(),
[f,epsilon](decltype(*dfreqs.cbegin()) x) {
return (std::abs(x-f) < epsilon);
})
;
if (found == dfreqs.cend()) { // no double frequencies
dfreqs.push_back (f);
expandFrequencies (f, n);
}
}
}
}

// no excitations
if (negfreqs.getSize () == 0) {
if (negfreqs.size () == 0) {
// use specified frequency
f = getPropertyDouble ("f");
dfreqs.add (f);
dfreqs.push_back (f);
expandFrequencies (f, n);
}

// build frequency dimension lengths
ndfreqs = new int[dfreqs.getSize ()];
for (i = 0; i < dfreqs.getSize (); i++) {
ndfreqs = new int[dfreqs.size ()];
for (i = 0; i < dfreqs.size (); i++) {
ndfreqs[i] = (n + 1) * 2;
}

Expand All @@ -419,17 +428,17 @@ void hbsolver::collectFrequencies (void) {
#endif /* HB_DEBUG */

// build list of positive frequencies including DC
for (n = 0; n < negfreqs.getSize (); n++) {
if ((f = negfreqs (n)) < 0.0) continue;
rfreqs.add (f);
for (n = 0; n < negfreqs.size (); n++) {
if ((f = negfreqs[n]) < 0.0) continue;
rfreqs.push_back (f);
}
lnfreqs = rfreqs.getSize ();
nlfreqs = negfreqs.getSize ();
lnfreqs = rfreqs.size ();
nlfreqs = negfreqs.size ();

// pre-calculate the j[O] vector
OM = new tvector<nr_complex_t> (nlfreqs);
for (n = i = 0; n < nlfreqs; n++, i++)
OM_(n) = nr_complex_t (0, 2 * pi * negfreqs (i));
OM_(n) = nr_complex_t (0, 2 * M_PI * negfreqs[i]);
}

// Split netlist into excitation, linear and non-linear part.
Expand Down Expand Up @@ -567,8 +576,8 @@ void hbsolver::createMatrixLinearA (void) {
A = new tmatrix<nr_complex_t> ((N + M) * lnfreqs);

// through each frequency
for (int i = 0; i < rfreqs.getSize (); i++) {
freq = rfreqs (i);
for (int i = 0; i < rfreqs.size (); i++) {
freq = rfreqs[i];
// calculate components' MNA matrix for the given frequency
for (auto *lc : lincircuits)
lc->calcHB (freq);
Expand Down Expand Up @@ -873,8 +882,8 @@ void hbsolver::calcConstantCurrent (void) {
circuit * vs = *it;
vs->initHB ();
vs->setVoltageSource (0);
for (int f = 0; f < rfreqs.getSize (); f++) { // for each frequency
nr_double_t freq = rfreqs (f);
for (int f = 0; f < rfreqs.size (); f++) { // for each frequency
nr_double_t freq = rfreqs[f];
vs->calcHB (freq);
VC (vsrc * lnfreqs + f) = vs->getE (VSRC_1);
}
Expand Down Expand Up @@ -920,7 +929,7 @@ int hbsolver::checkBalance (void) {
nr_double_t iabstol = getPropertyDouble ("iabstol");
nr_double_t vabstol = getPropertyDouble ("vabstol");
nr_double_t reltol = getPropertyDouble ("reltol");
int n, len = FV->getSize ();
int n, len = FV->size ();
for (n = 0; n < len; n++) {
// check iteration voltages
nr_double_t v_abs = abs (VS->get (n) - VP->get (n));
Expand Down Expand Up @@ -1081,13 +1090,15 @@ void hbsolver::loadMatrices (void) {
}

/* The following function transforms a vector using a Fast Fourier
Transformation from the time domain to the frequency domain. */
Transformation from the time domain to the frequency domain.
\todo rewrite ugly sould die
*/
void hbsolver::VectorFFT (tvector<nr_complex_t> * V, int isign) {
int i, k, r;
int n = nlfreqs;
int nd = dfreqs.getSize ();
int nodes = V->getSize () / n;
nr_double_t * d = (nr_double_t *) V->getData ();
int nd = dfreqs.size ();
int nodes = V->size () / n;
nr_double_t * d = (double *)V->getData ();

if (nd == 1) {
// for each node a single 1d-FFT
Expand Down Expand Up @@ -1316,7 +1327,7 @@ void hbsolver::fillMatrixLinearExtended (tmatrix<nr_complex_t> * Y,
int nnode = vs->getNode(NODE_2)->getNode ();
for (int f = 0; f < lnfreqs; f++, sc++) { // for each frequency
// fill right hand side vector
nr_double_t freq = rfreqs (f);
nr_double_t freq = rfreqs[f];
vs->calcHB (freq);
I_(sc) = vs->getE (VSRC_1);
// fill MNA entries
Expand Down Expand Up @@ -1390,7 +1401,7 @@ void hbsolver::saveResults (void) {
}
// save frequency vector
if (runs == 1) {
for (int i = 0; i < lnfreqs; i++) f->add (rfreqs (i));
for (int i = 0; i < lnfreqs; i++) f->add (rfreqs[i]);
}
// save node voltage vectors
int nanode = 0;
Expand Down
10 changes: 6 additions & 4 deletions qucs-core/src/hbsolver.h
Original file line number Diff line number Diff line change
Expand Up @@ -25,6 +25,8 @@
#ifndef __HBSOLVER_H__
#define __HBSOLVER_H__

#include <vector>

#include "ptrlist.h"
#include "tvector.h"

Expand Down Expand Up @@ -85,11 +87,11 @@ class hbsolver : public analysis
void saveNodeVoltages (circuit *, int);

private:
tvector<nr_double_t> negfreqs; // full frequency set
tvector<nr_double_t> posfreqs; // full frequency set but positive
tvector<nr_double_t> rfreqs; // real positive frequency set
std::vector<nr_double_t> negfreqs; // full frequency set
std::vector<nr_double_t> posfreqs; // full frequency set but positive
std::vector<nr_double_t> rfreqs; // real positive frequency set
int * ndfreqs; // number of frequencies for each dimension
tvector<nr_double_t> dfreqs; // base frequencies for each dimension
std::vector<nr_double_t> dfreqs; // base frequencies for each dimension
nr_double_t frequency;
strlist * nlnodes, * lnnodes, * banodes, * nanodes, * exnodes;
ptrlist<circuit> excitations;
Expand Down
6 changes: 3 additions & 3 deletions qucs-core/src/interface/e_trsolver.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -371,7 +371,7 @@ void e_trsolver::printx()
{
char buf [1024];

for (int r = 0; r < x->getSize(); r++) {
for (int r = 0; r < x->size(); r++) {
buf[0] = '\0';
//sprintf (buf, "%+.2e%+.2ei", (double) real (x->get (r)), (double) imag (x->get (r)));

Expand Down Expand Up @@ -828,9 +828,9 @@ void e_trsolver::copySolution (tvector<nr_double_t> * src[8], tvector<nr_double_
for (int i = 0; i < 8; i++)
{
// check sizes are the same
assert (src[i]->getSize () == dest[i]->getSize ());
assert (src[i]->size () == dest[i]->size ());
// copy over the data values
for (int j = 0; j < src[i]->getSize (); j++)
for (int j = 0; j < src[i]->size (); j++)
{
dest[i]->set (j, src[i]->get (j));
}
Expand Down
2 changes: 1 addition & 1 deletion qucs-core/src/nasolver.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -248,7 +248,7 @@ void nasolver<nr_type_t>::applyNodeset (bool nokeep)
if (x == NULL || nlist == NULL) return;

// set each solution to zero
if (nokeep) for (int i = 0; i < x->getSize (); i++) x->set (i, 0);
if (nokeep) for (int i = 0; i < x->size (); i++) x->set (i, 0);

// then apply the nodeset itself
for (nodeset * n = subnet->getNodeset (); n; n = n->getNext ())
Expand Down
57 changes: 41 additions & 16 deletions qucs-core/src/spline.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -30,6 +30,7 @@
#include <stdlib.h>
#include <string.h>
#include <assert.h>
#include <vector>

#include "logging.h"
#include "complex.h"
Expand Down Expand Up @@ -59,7 +60,17 @@ spline::spline (int b) {
}

// Constructor creates an instance of the spline class with vector data given.
spline::spline (vector y, vector t) {
spline::spline (qucs::vector y, qucs::vector t) {
x = f0 = f1 = f2 = f3 = NULL;
d0 = dn = 0;
n = 0;
boundary = SPLINE_BC_NATURAL;
vectors (y, t);
construct ();
}

// Constructor creates an instance of the spline class with tvector data given.
spline::spline (::std::vector<nr_double_t> y, ::std::vector<nr_double_t> t) {
x = f0 = f1 = f2 = f3 = NULL;
d0 = dn = 0;
n = 0;
Expand All @@ -82,7 +93,7 @@ spline::spline (tvector<nr_double_t> y, tvector<nr_double_t> t) {
#define y_ (y)

// Pass interpolation datapoints as vectors.
void spline::vectors (vector y, vector t) {
void spline::vectors (qucs::vector y, qucs::vector t) {
int i = t.getSize ();
assert (y.getSize () == i && i >= 3);

Expand All @@ -93,10 +104,22 @@ void spline::vectors (vector y, vector t) {
}
}

// Pass interpolation datapoints as tvectors.
void spline::vectors (::std::vector<nr_double_t> y, ::std::vector<nr_double_t> t) {
int i = (int)t.size ();
assert ((int)y.size () == i && i >= 3);

// create local copy of f(x)
realloc (i);
for (i = 0; i <= n; i++) {
f0[i] = y[i]; x[i] = t[i];
}
}

// Pass interpolation datapoints as tvectors.
void spline::vectors (tvector<nr_double_t> y, tvector<nr_double_t> t) {
int i = t.getSize ();
assert (y.getSize () == i && i >= 3);
int i = t.size ();
assert (y.size () == i && i >= 3);

// create local copy of f(x)
realloc (i);
Expand Down Expand Up @@ -213,7 +236,7 @@ void spline::construct (void) {
// second kind of cubic splines
else if (boundary == SPLINE_BC_PERIODIC) {
// non-trigdiagonal equations - periodic boundary condition
nr_double_t * z = new nr_double_t[n+1];
::std::vector<nr_double_t> z (n+1);
if (n == 2) {
nr_double_t B = h[0] + h[1];
nr_double_t A = 2 * B;
Expand All @@ -227,18 +250,20 @@ void spline::construct (void) {
}
else {
tridiag<nr_double_t> sys;
tvector<nr_double_t> o (n);
tvector<nr_double_t> d (n);
tvector<nr_double_t> b;
b.setData (&z[1], n);
::std::vector<nr_double_t> o (n);
::std::vector<nr_double_t> d (n);
::std::vector<nr_double_t> b(&z[1],&z[n]);
//b.setData (&z[1], n);
for (i = 0; i < n - 1; i++) {
o(i) = h[i+1];
d(i) = 2 * (h[i+1] + h[i]);
b(i) = 3 * ((f0[i+2] - f0[i+1]) / h[i+1] - (f0[i+1] - f0[i]) / h[i]);
o[i] = h[i+1];
d[i] = 2 * (h[i+1] + h[i]);
b[i] = 3 * ((f0[i+2] - f0[i+1]) / h[i+1] - (f0[i+1] - f0[i]) / h[i]);
z[i+1] = b[i];
}
o(i) = h[0];
d(i) = 2 * (h[0] + h[i]);
b(i) = 3 * ((f0[1] - f0[i+1]) / h[0] - (f0[i+1] - f0[i]) / h[i]);
o[i] = h[0];
d[i] = 2 * (h[0] + h[i]);
b[i] = 3 * ((f0[1] - f0[i+1]) / h[0] - (f0[i+1] - f0[i]) / h[i]);
z[i+1] = b[i];
sys.setDiagonal (&d);
sys.setOffDiagonal (&o);
sys.setRHS (&b);
Expand All @@ -248,7 +273,7 @@ void spline::construct (void) {
}

f1 = new nr_double_t[n+1];
f2 = z; // reuse storage
f2 = &z.front (); // reuse storage
f3 = h;
for (i = n - 1; i >= 0; i--) {
f1[i] = (f0[i+1] - f0[i]) / h[i] - h[i] * (z[i+1] + 2 * z[i]) / 3;
Expand Down
5 changes: 4 additions & 1 deletion qucs-core/src/spline.h
Original file line number Diff line number Diff line change
Expand Up @@ -25,6 +25,7 @@
#ifndef __SPLINE_H__
#define __SPLINE_H__

#include <vector>
#include "tvector.h"

namespace qucs {
Expand All @@ -45,12 +46,14 @@ class spline
public:
spline ();
spline (int);
spline (qucs::vector, qucs::vector);
spline (tvector<nr_double_t>, tvector<nr_double_t>);
spline (qucs::vector, qucs::vector);
spline (::std::vector<nr_double_t>, ::std::vector<nr_double_t>);
~spline ();

void vectors (qucs::vector, qucs::vector);
void vectors (tvector<nr_double_t>, tvector<nr_double_t>);
void vectors (::std::vector<nr_double_t>, ::std::vector<nr_double_t>);
void vectors (nr_double_t *, nr_double_t *, int);
void construct (void);
poly evaluate (nr_double_t);
Expand Down
4 changes: 2 additions & 2 deletions qucs-core/src/tmatrix.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -289,7 +289,7 @@ tmatrix<nr_type_t> operator * (tmatrix<nr_type_t> a, tmatrix<nr_type_t> b) {
// Multiplication of matrix and vector.
template <class nr_type_t>
tvector<nr_type_t> operator * (tmatrix<nr_type_t> a, tvector<nr_type_t> b) {
assert (a.getCols () == b.getSize ());
assert (a.getCols () == b.size ());
int r, c, n = a.getCols ();
nr_type_t z;
tvector<nr_type_t> res (n);
Expand All @@ -304,7 +304,7 @@ tvector<nr_type_t> operator * (tmatrix<nr_type_t> a, tvector<nr_type_t> b) {
// Multiplication of vector (transposed) and matrix.
template <class nr_type_t>
tvector<nr_type_t> operator * (tvector<nr_type_t> a, tmatrix<nr_type_t> b) {
assert (a.getSize () == b.getRows ());
assert (a.size () == b.getRows ());
int r, c, n = b.getRows ();
nr_type_t z;
tvector<nr_type_t> res (n);
Expand Down

0 comments on commit d983fc5

Please sign in to comment.