• Welcome to the new Internet Infidels Discussion Board, formerly Talk Freethought.

The Programming Thread

Bad bet. APL is much more efficient. My little brother, almost 78, up in Ilwaco WN has been replacing COBOL with APL for banks, Insurance companies, and Oil Companies for over 40 years.
APL? There's a programming language with that name, but it uses lots of extra symbols, and it is very terse. Why APL and not (say) Python?

Because Python was not around when New Jersey digitized their unemployment claim processing in the 1970s. They had a choice between COBOL, Lisp, ForTran, and APL, and a few more that are even more obscure today - and two of these explicitly targeted a scientific rather than administrative applications.
 
I’m not doing anything serious, but needed to plot for loking at math and sciene. I didn’t want to go to the trouble of installing a graphics package so I took the path of least resistance and imported to Scilab. Scilab is ok for coding but it is a slow intrpreter and the debugger is not very good.
Visual Studio is pretty good.

Sciab has a full range of 2d.3d plots abd exort functions. Easy install.

#include "steve.h"
#include <fstream>
#include <iomanip>
using namespace std;

int main()
{
ofstream testfile;
testfile.open("example.txt");
float x,y, dt,t;
dt = 0.01;
for (t = 0; t <= 1 ; t += dt) {
x = sin(_PI2 * t);
y = cos(_PI2 * t);
cout << t << setw(20) << x << setw(20) << y << "\n";
testfile << t << setw(20) << x << setw(20) << y << "\n";
}
testfile.close();
return 0;
}
//sclab
clear
[f1,err] = mopen("c:\\documents\SCILAB\example.txt","rt")
if(err < 0)then disp("file open error");end;
// -1 go to end of file
x=mfscanf(-1,f1,' %e %e %e')
// r rows c columns
[r,c] = size(x)
err = mclose("all")
if(err < 0) then disp("file close error")end;
plot2d(x(1:r,1),x(1:r,2))
g=gca()
g.x_location = "middle"
xs2png(w1,"c:\documents\scilab\plots\fig1")

fig1.png
 
A make work project to keep busy, a low pass digital filter. It us under the category of real time systems.

In a real implementation an analog to digital converter would be sampling a signal and pushing samples into a buffer the length of the number of filter coefficients. The filter has to execute and an output generated before the next input sample.

The filter is auto-regressive. At any sample n the filter operates on the sample plus the previous number of samples equal to the filter length.

There are different algorithms to generate filter coefficients. That would be a math thread.

The simulation is mot entirely accurate, it uses floating point. In a real system the output and input to analog digital converters are signed integers that limit the accuracy and resolution.



#include "steve.h""
#include <fstream>
#include <iomanip>
using namespace std;

void make_sin(float fs, double x[], double y[],int sf, double tmax, int msize) {
int i;
double dt,t;
dt = 1/fs;
t = 0;
for (i = 0; i < msize; i++) {
x = t;
y = sin(_PI2 * sf * t);
t += dt;
}
}

void fir_lp(double isig[], double osig[], double coef[], int lencoef, int lensig) {
// low pass filter
int i, j, k, n;
double acc;
double* buffer = new double[lencoef];
// initialize the input buffer
for (i = 0; i < lencoef; i++)buffer = 0;
k = 0;
for (i = 0; i < lensig; i++) {
acc = 0;
for (j = 0; j < lencoef; j++) acc += coef[j] * buffer[j]);
osig[k++] = acc; // output to digital to analog converter
//push new analog to digital converter sample onto the buffer.
for (n = 1; n < lencoef; n++)buffer[n - 1] = buffer[n];
buffer[lencoef - 1] = isig;
}
delete [] buffer;
}

int main()
{
double samp_freq,tmax;
int sig_freq, i, n, ncoef, msize;
// get filter coefficients
ifstream coeffile;
coeffile.open("coef.txt");
coeffile >> ncoef; // filter length
coeffile >> tmax; // simulation time
double* coef = new double[ncoef];
for (i = 0; i < ncoef; i++) coeffile >> coef;
coeffile.close();
sig_freq = 1; //hertz sampling frequency
samp_freq = 100; //hertz sample frequency
msize = ceil(tmax * samp_freq); //array size
double* x = new double[msize];
double* filtin = new double[msize];
double* filtout = new double[msize];
ofstream outfile;
outfile.open("lpout.txt");
make_sin(samp_freq, x, filtin, sig_freq,tmax,msize);
fir_lp(filtin, filtout, coef,ncoef,msize);
for(i=0;i < msize; i++)outfile << x << " " << filtin << " "<< filtout << "\n";
outfile.close();
delete [] x;
delete [] filtin;
delete [] filtout;
return(0);
}
 
....
Code:
#include "steve.h""
#include <fstream>
#include <iomanip>
using namespace std;

void make_sin(float fs, double x[], double y[],int sf, double tmax, int msize) {
    int i;
    double dt,t;
    dt = 1/fs;
    t = 0;
    for (i = 0; i < msize; i++) {
        x[i] = t;
        y[i] = sin(_PI2 * sf * t);
        t += dt;
    }
}
Hi I was just testing out the code tag on your code.... it even automatically fixed all of the indentation.... which is nice....

Code:
[COLOR=#333333]int main()[/COLOR]
[COLOR=#333333]{[/COLOR]
[COLOR=#333333]double samp_freq,tmax;[/COLOR]
[COLOR=#333333]int sig_freq, i, n, ncoef, msize;[/COLOR]
[COLOR=#333333]// get filter coefficients[/COLOR]
[COLOR=#333333]ifstream coeffile;[/COLOR]
[COLOR=#333333]coeffile.open("coef.txt");[/COLOR]
[COLOR=#333333]coeffile >> ncoef; // filter length[/COLOR]
[COLOR=#333333]coeffile >> tmax; // simulation time[/COLOR]
[COLOR=#333333]double* coef = new double[ncoef];[/COLOR]
[COLOR=#333333]for (i = 0; i < ncoef; i++) coeffile >> coef[i];[/COLOR]
[COLOR=#333333]coeffile.close();[/COLOR]
[COLOR=#333333]sig_freq = 1; //hertz sampling frequency[/COLOR]
[COLOR=#333333]samp_freq = 100; //hertz sample frequency[/COLOR]
[COLOR=#333333]msize = ceil(tmax * samp_freq); //array size[/COLOR]
[COLOR=#333333]double* x = new double[msize];[/COLOR]
[COLOR=#333333]double* filtin = new double[msize];[/COLOR]
[COLOR=#333333]double* filtout = new double[msize];[/COLOR]
[COLOR=#333333]ofstream outfile;[/COLOR]
[COLOR=#333333]outfile.open("lpout.txt");[/COLOR]
[COLOR=#333333]make_sin(samp_freq, x, filtin, sig_freq,tmax,msize);[/COLOR]
[COLOR=#333333]fir_lp(filtin, filtout, coef,ncoef,msize);[/COLOR]
[COLOR=#333333]for(i=0;i < msize; i++)outfile << x[i] << " " << filtin[i] << " "<< filtout[i] << "\n";[/COLOR]
[COLOR=#333333]outfile.close();[/COLOR]
[COLOR=#333333]delete [] x;[/COLOR]
[COLOR=#333333]delete [] filtin;[/COLOR]
[COLOR=#333333]delete [] filtout;[/COLOR]
[COLOR=#333333]return(0);[/COLOR]
[COLOR=#333333]}[/COLOR]
Edit: it didn't fix the indentation for that - oh well....
 
The frequency response of the filter and the output for a sin input. The startup delay is the buffer filling up.


firmag.png

firout.png
 
The coef .txt file

33
2
-0.0086827058
0.0144278375
0.0113203161
-0.0080360033
-0.0091899954
-0.0009891973
0.0227447187
0.0087521368
-0.0368257634
-0.0127543525
0.0378075284
0.0365069835
-0.0368261824
-0.0971619203
0.0483744099
0.3187605599
0.4425552553
0.3187605599
0.0483744099
-0.0971619203
-0.0368261824
0.0365069835
0.0378075284
-0.0127543525
-0.0368257634
0.0087521368
0.0227447187
-0.0009891973
-0.0091899954
-0.0080360033
0.0113203161
0.0144278375
-0.0086827058
 
Found this on the net while looking for something else. Truncation/rounding at a decimal place. Shift left, round to get rid of unwanted digits, shift back.


double dec_places(double val, unsigned long long int n) {
unsigned long long int dp = pow(10, n);
return(round(val * dp) / dp);
}

int main() {
double xf, x = 1.123456789;
unsigned long long int n_dec = 5;
xf = dec_places(x, n_dec);
printf(" x %2.15f xf %2.15f \n", x, xf);
}
 
That leads to think. What might be the best sorts of languages for exploratory / experimental sorts of programming?

I like computer-algebra software.

 List of computer algebra systems
 
If you want to do number crunching with big blocks of data, here is what might be useful:

For general-purpose programming, I propose this spectrum:
  • Python - dynamic typing, easy to start with
  • Java - static typing, generics, garbage collection
  • C++ - static typing, generics, semi-explicit memory management - designed for high performance
For experimental or exploratory programming, I do NOT recommend plain C. One has to do a LOT of stuff by hand by it, and its preprocessor and its memory management are just plain bug-prone.

C++ supports several safe alternatives.

One can declare constants with "const", as a type-safe alternative to preprocessor declarations.

One can also use "const" for function args and method declarations, as a way of stating that something will not be changed.

Templates are a good type-safe alternative to preprocessor functions.

Its Standard Template Library is a set of collection classes and functions for working with them, and a nice thing about them is that one does not have to worry about memory leaks. The collections deallocate all their allocated memory when they go out of scope.

C++ now supports anonymous functions (lambdas), and that should make some uses of it a LOT easier.
 
I used Mathcad for a long time. It has a white board GUI. Drag and place symbols and it figures it out. Popu;ar in engnerering.

I prefer typed languages and C.Packages like Scilab do everything dynamically. Define an array or variable by using it. Its great for quick analysis to a point. For me pointers in C were essential when working with multiple data sets.

If oiu do a lot of work with matrices having matrix functions is a help.

Matlab has been a standard since the 80s in engineering. Toolboxes like control systems. It is known for simulation. Matlab is a good tool especially in organization where code can be shared and nobody wants to be a programmer.

In C/C++ you can create numbers of any length limited by memory.

Scilab can invoke an installed C compiler and compile Scilab code to a library, make your own functions. You can create sophisticated apps with windows and pull down menus.I am no longer doing actual work. I use VS C++ and export to Scilab fpr graphics

If I were to use a simulation tool at work it would be Matlab if the company paid for it . Scilab is a competitor and it is free. I carried it around from job to job.
 
Object Oriented Programming.

Calculus is not limited to physical sciences, derivatives and related rates are used in business and economics.

I grew to approach code like hardware. Break the problem down into manageable parts.

Come up with a set of objects to provide functionality on which to build the end application. Bild and test an object one function at tme testing each function and then the finished object. Once the object is working it can be spawned to accommodate a number of data sets.

The object takes data which can be from anything and calculates the first two derivatives. The data and table are protected from code outside the object. The only access is through a function and public variables.


class difference_table {
public:
string data_file, output_file;
int error = 0;
double y_out = 0, dy1_out = 0, dy2_out = 0, x_in = 0;
void create_table(void);
void load_data(void);
void get_value(void);
void save(void);

private:
int n;
double* x = new double[2000];
double* y = new double[2000];
double* dy1 = new double[2000];
double* dy2 = new double[2000];
};

void difference_table::get_value(void) {
//query the table, if no exact match linear interpolation
int i;
double k;
// check if query within data range
if (x_in < x[0] || x_in > x[n-1]) { error = 1; return; }
for(i = 1;i<n-1;i++){
// check for exact match
if (x_in == x) {
y_out = y;
dy1_out = dy1;
dy2_out = dy2;
return;
}
// interpolate
if((x_in > x) && (x_in < x[i+1]) && ( i < n-2)){
k = (x_in - x) /(x[i + 1] - x);
y_out = y + ((y[i + 1] - y) * k);
dy1_out = dy1 + ((dy1[i + 1] - dy1) * k);
dy2_out = dy2 + ((dy2[i + 1] - dy2) * k);
}
}
}//get_value()

void difference_table::save(void) {
int i;
ofstream datafile;
datafile.open(output_file);
datafile << n << "\n";
for (i = 0; i < n; i++)
datafile << x << " " << y << " " << dy1 << " " << dy2 << "\n";
datafile.close();
}//save()


void difference_table::create_table(void) {
// dy/dx 1st and dy'/dx 2nd dervitve of the data
int i;
double dx;
for (i = 0; i < n - 2; i++){
dx = x[i + 1] - x;
dy1 = (y[i + 1] - y) / dx;
}
for (i = 0; i < n - 3; i++) {
dx = x[i + 1] - x;
dy2 = (dy1[i + 1] - dy1) / dx;
}
}//create_table()


void difference_table::load_data(void){
int i;
for (i = 0; i < 2000; i++) { y = 0; dy1 = 0; dy2 = 0; }
ifstream datafile;
datafile.open(data_file);
datafile >> n ;
for (i = 0; i < n; i++)datafile >> x >> y;
datafile.close();
}//load_data()
 
Testing with a sin wave

void main() {

double* y = new double[2000];
double* x = new double[2000];
double t = 0, tmax = 1, dt = .001;
int i, n;

difference_table t1;
//sin wave test data
n = floor(tmax / dt);
for (i = 0; i < n; i++){
y = sin(_PI2 * t);
x = t;
t += dt;
}
save("data.txt", n - 1, x, y);

t1.data_file = "data.txt";
t1.output_file = "table_out.txt";
t1.load_data();
t1.create_table();
t1.save();
t1.x_in = .0625;
t1.get_value();
cout << t1.y_out;

}
sin in
yin.png
1st derivative
dy1.png
2nd derivative
dy2.png
 
Ever wonder how numbers are multiplied in a PC?


It is done at thebit level.


An nxn binarymultipy requires a 64 bit accumulator. The algorithm is as withdecimal multiplication. Multiply, shift, and add. In hardware it iscalled a MAC, multiplier accumaulatr. A barrel shifter is a fastasynchrounous, non clocked, shift register.


In binary it iseasy. If the mutiplier bit is 1 add the multiplcand to anaccumulator. Shift and repeat.


When multipyingfrcational numbers the rrsult in the accumulatir has to be shifted nbits to te right, analogus to setting the decimal point in decimalmultiplication.


srart number =.0095810000 – compiler mult
16 bit mult =.0092851593
32 bit mult =.0095809999


unsigned long longmult(unsigned long long a,unsigned long long b,int nbits) {
// 32x32 maxunsigned binary multiply
int i;
unsigned long long mask = 1, acc = 0;
for (i = 0; i <nbits; i++)if (a & (mask << i)) acc += (b << i);
return(acc);
}//mult()


unsigned long long dec2bin( int nbits, double num) {
// finds a binary approx to a num < 1 based on nbits
double lsb = 1 /pow(2, nbits);
return(ceil(num /lsb));
}//dec2bin()


doublefrac2float(unsigned long long ,int nbits) {
//cnvertsfractiianl binary to decimal
double lsb = 1 /pow(2,nbits);
return(lsb * a);
}//frac2float()


void main(void) {
unsigned long longint i1, i2, nb, int_result;
double num_start,fl1, fl2, dec_result, lsb;
nb = 16;
// fractional mult
fl1 = .13; //floats to multiply
fl2 = .0737;
num_start = fl1 *fl2;
// convert tointeger
i1 = dec2bin(nb,fl1);
i2 = dec2bin(nb,fl2);
//frac mult
int_result = mult(i1,i2,nb);
int_result=int_result >> nb;
dec_result =frac2float(int_result,nb);
printf("num_start%1.11f dec mult a*b %1.11f\n", num_start, dec_result);
}//main()
 
https://www.geeksforgeeks.org/introduction-of-k-map-karnaugh-map/
https://en.wikipedia.org/wiki/Karnaugh_map




Bit integer addition. An exercise in formal logic.
a,b integer inputs
s integer sum
c carry in
cout carry out
c = 0
0 + 0 + c = 0
0 + 1 + c = 1
1 + 0 + c = 1
1 + 1 + c = 0 carry1
c = 1
0 + 0 + c = 1
0 +1 + c = 0 carry1
1 + 0 + c = 0 carry1
1 + 1 + c = 1 carry1


Truth table
a b c s cout
0 0 0 0 0
0 1 0 1 0
1 0 0 1 0
1 1 0 0 1


0 0 1 1 0
0 1 1 0 1
1 0 1 0 1
1 1 1 1 1


K Maps are usefulmanual tool to dervive lagical expwerssons for comlicated logicfunctions. A truth table can represent anything, such as makingdecsions based on data base items. T and F can be used instead of 1and 0;


A brute forceappoach might be to write a logic expresson for each line in thetable and OR them.


Without resorting toK Maps by inspection.


S is true when a orb is true but not both without carry, XOR Exclusive OR .With carry it is XOR inverted.


s = (a^b) &&!c) || (!(a^b) && c))

cout is true for aand b true regardles of carry. The other cout true cinditions are (aXOR b ) AND c

cout = (a&b) ||((a^b) && c)

class bmath {
public:
int n_bits;
double dec_result,num;
unsigned long longresult, a, b;


void mult(void);
voidfrac2float(void);
void dec2bin(void);
void sub(void);
void badd(void);
};


void bmath::badd(void) {
int i;
unsigned long longmask = 1, carry = 0, c = 0, bita = 0, bitb = 0;
for (i = 0; i <n_bits; i++) {
bita = a &mask;
bitb = b &mask;
if (((bita ^ bitb)&& !carry) || (!(bita ^ bitb) && carry))c = c | mask;
if ((bita &&bitb) || ((bita ^ bitb) && carry)) carry = 1; else carry = 0;
mask <<= 1;
}//for
result = c;
}//badd()


void main(void) {
unsigned long long,a,b,x;
bmath m1;
m1.n_bits = 32
m1.a = 7;
m1.b = 7;
m1.badd();
printf(" a %lld b %lld trdult = %ll x = %lld \n", m1.a, m1.b,m1.result,x);
}//main()
 
Back in the 20th-century I was a programmer with almost a fetish for time-efficient low-level code. Multiplication and division were often important bottlenecks. I'll content myself with just the most interesting example of what can be done:
Code:
#include        <stdio.h>
#include        <stdlib.h>

/* Use this approach when there is a smallish known bound on multiplier size */
#define BOUND   20000

unsigned int sq[BOUND * 4];
#define osq (sq + BOUND * 2)

/* This multiplier is FASTER than '*' on almost all processors without special assist */
int mult(int a, int b)
{
        return osq[a + b] - osq[a - b];
}

int main(int argc, char **argv)
{
        int a, b;

        /* Setup sq table */
        for (a = 0; a < BOUND * 2; a++)
                osq[-a] = osq[a] = a * a / 4;

        /* Verify mult() works by testing every case */
        for (a = 0; a < BOUND; a++)
        for (b = 0; b < BOUND; b++)
                if (a * b != mult(a, b))
                        printf("Mult failed at %d %d\n", a, b);
        exit(0);
}
Do you see the cute algebraic trick at work here? I have verified that when it is applicable, this code often outperforms a 'normal' multiply. (The example code assumes a,b non-negative, but can be adapted to support negative numbers also.)


BTW, the Intel 80x86 has a single-instruction multiply by 5: leal (%eax,%eax,4), %eax. The Motorola 680xx has a similar instruction.
 
Back in the 20th-century I was a programmer with almost a fetish for time-efficient low-level code. Multiplication and division were often important bottlenecks. I'll content myself with just the most interesting example of what can be done:
Code:
#include        <stdio.h>
#include        <stdlib.h>

/* Use this approach when there is a smallish known bound on multiplier size */
#define BOUND   20000

unsigned int sq[BOUND * 4];
#define osq (sq + BOUND * 2)

/* This multiplier is FASTER than '*' on almost all processors without special assist */
int mult(int a, int b)
{
        return osq[a + b] - osq[a - b];
}

int main(int argc, char **argv)
{
        int a, b;

        /* Setup sq table */
        for (a = 0; a < BOUND * 2; a++)
                osq[-a] = osq[a] = a * a / 4;

        /* Verify mult() works by testing every case */
        for (a = 0; a < BOUND; a++)
        for (b = 0; b < BOUND; b++)
                if (a * b != mult(a, b))
                        printf("Mult failed at %d %d\n", a, b);
        exit(0);
}
Do you see the cute algebraic trick at work here? I have verified that when it is applicable, this code often outperforms a 'normal' multiply. (The example code assumes a,b non-negative, but can be adapted to support negative numbers also.)


BTW, the Intel 80x86 has a single-instruction multiply by 5: leal (%eax,%eax,4), %eax. The Motorola 680xx has a similar instruction.

It is just a freaking exercise.

All processors had/have signed and unsigned integer addition ADD ADDC, and hardware multiply. On the 8086 motherboards there was a socket for a hardware math chip that operated on floating point numbers directly. It was integrated on one of the x86 processors into the chip.

To get larger numbers you had to chain 8 bit bytes using an add with carry instruction. A 16 bit addition would take to 8 bit adds.

It is no trick. It is a routine and and simple exercise in logic minimization.

When deigning hardware I used truth tables. Beyond a certain level it is impossible to do in the head. K maps and other logic minimization techniques are now automated. Tools take inputs and generate logic.

Boolean Algebra is a a form of formal logic and is routine in engineering.

BTW, since the 80s nobody designs IC scale logic using direct digital logic, it is done in software like VERILOG that takes a c like language and creates logic to emulate the code.
 
... I'll content myself with just the most interesting example of what can be done:
...

It is just a freaking exercise.

All processors had/have signed and unsigned integer addition ADD ADDC, and hardware multiply.

I think you're correct that essentially all of today's processors have very fast multiplication. I don't think that detracts from possible interest in the code I posted.

And I think your 'had/have' should be amended. As recently as 28 years ago (when last I indulged in the most extreme forms of my efficiency fetish!), the code I posted above multiplied MUCH faster than the default multiplication on a Sun SPARCstation, even when the assembly-code for .mul was replaced with a faster version I wrote.
 
You are oversimplifying .

As I said before, it comes dwn to the assembly code generated by the compiler, and not all compilers gnerate the exact same assembly code. Fastest math code is always done in assembly code.

As I said before, you would have to use a debugger/simulations to determine number of clock cycles for different methods. That is the 'proof'.

Multiple loops require multiple reloading and initialing of the loop counter in the chip. And so on. Multiple fetches to memory.

These days considering all the man hours gone into math optimization I doubt there is much in the way of improvements.

What incresed speed was obviosly clock speed, but also memory width. Wider memory meant fewer fetches from memory' On an 8 bit machine two fetches were needed for a 16 bit integer. On a 64 bit machine a 64 bit int takes one fetch.

void main(void) {
unsigned long long ,a,b,x;
bmath m1;
m1.n_bits = 32
m1.a = 7;
m1.b = 7;
m1.badd();
printf(" a %lld b %lld trdult = %ll x = %lld \n", m1.a, m1.b, m1.result,x);
}//main()

This was my secomd version

void bmath::badd(void) {
unsigned int nbits = 32,i, carry = 0;
unsigned long long mask = 1, c = 0,bita,bitb;

for (i = 0; i < n_bits; i++) {
bita = mask & a;
bitb = mask & b;
// 1+0,0+1
if(bita ^ bitb){if(!carry)c = c|mask;else carry = 1;}
//0+0
else if(!bita & !bitb){if(carry)c = c|mask; carry = 0;}
//1+1
else if(bita && bitb){if(carry)c = c|mask; carry = 1;}
mask <<= 1;
}//for
result = c;
}//badd()

My first version was a brute force line by line implantation of the truth table.
 
You are oversimplifying .

As I said before, it comes dwn to the assembly code generated by the compiler, and not all compilers gnerate the exact same assembly code. Fastest math code is always done in assembly code.

As I said before, you would have to use a debugger/simulations to determine number of clock cycles for different methods. That is the 'proof'.

You do not have the slightest idea whom you are talking to! Do you think I make such claims without testing speed and examining compiled code? :confused:

There's a good chance I was disassembling and optimizing code when you were still in diapers. Have I mentioned that at least two Engineering Vice Presidents once called me "the best microprogrammer in Silicon Valley"?
 
Back
Top Bottom