Show Posts
Pages: 1 ... 7 8 [9] 10 11 ... 25
121  Using Arduino / General Electronics / Re: Getting impulses (0.5 Wh) from an home energy counter on: December 01, 2013, 06:12:23 pm
It appears that this device is produced by a Dutch company called "inepro."  The manufacturer intends for resellers to apply their own brand names to this product.  Here's a link to a French version of what seems to be the datasheet:  http://www.bis-electric.com/catalog/doc_technique/M32M/NOTICE_M32M.pdf. I think that the pulse output, called, "Sortie d'impulsions," is on page 9.  An online translation program says that the pulse output is electrically isolated from the rest of the device, that the external source voltage must not exceed 27V, and that the maximum current is 27 mA.

Here's a link to an English version:  www.bis-electric.com/catalog/images/manuel_m32L.pdf.  It says the same thing, in English.  
122  Community / Exhibition / Gallery / Re: Arduino generating prime numbers on: November 26, 2013, 01:13:55 am
Eleven bad.  That should be a one.  If that statement ever executed, it should have caused an error; nevertheless, a test with five million primes ran without a hitch last night.

Here's what I did to test it:
  • Added a "print something and stop" inside those brackets.  The Due went to 120K primes found, and that didn't execute.
  • Generated a whole bunch of uint64_t random pairs, multiplied them, and checked about 40K results in Perl with Math::GMPz.  Again, the print and stop didn't execute.  All the products were correct.
  • Ran 100K random pairs, this time with the most significant nybble set to 0xF.  Again, no print and stop, and products checked out.
  • Rewrote the routine to operate on bytes, work on nybbles, and generate integer products; ran it against all 65K combinations of inputs, and checked.  Again, no print and stop, and all the products checked.
  • Multiplied 0xFFFFFFFFFFFFFFFF by itself.  Nothing went wrong.
Apparently, that statement never executes - there's never a carry after the second addition to the intermediate product.  That product captures the products of the high half of one input multiplied by the low half of the other.  That happens twice, and those products are added, and that may generate a carry - that's the first test.  Then the  high half of the product of the low halves of the inputs is added to the result.  I'd expect that addition to sometimes generate a carry, but it appears that it never does.  I can see why it won't when the inputs are restricted to 54 bits, but I can't see how it never gets a carry with arbitrary inputs.  If I can figure out why that test is unnecessary, I can remove it.

[Edit: Changed "11ULL" to "1ULL" in file attached to reply #87; fixed spelling.]
123  Community / Exhibition / Gallery / Re: Arduino generating prime numbers on: November 25, 2013, 04:35:19 pm
Indeed, the Miller-Rabin test calculates a whole lot of (x*y) % n's, where n is the number under test.  There's benefit to be had from doing something first with n if it will simplify the mod operations.

Here's "Barrett reduction," and it does just that, though without making a table: http://en.wikipedia.org/wiki/Barrett_reduction.  It lets us calculate an integer function of the number under test, about the same size, and use it to get modulo results by multiplication, shifting, and subtraction.  It adds complexity, but reduces the effort of multiple modulo calculations with the same modulus.  The problem is that for this project, the multiplication yields a 3-uint64_t product. Routines for multiple uint64_t's are inefficient compared to calculating with native data types, as they can't see the carry bit.  Barrett provides a speed improvement, but not the whopping acceleration that I hoped for.  It takes about 60% as long to calculate a 16-digit prime, average, as  the version from reply #81.  It beats one second per prime, average, on the Uno.

Barrett numbers k and m are added the cheesy way - as global variables, just before the Miller-Rabin tests.  There are a bunch of tweaks that would improve the speed a little bit, but the ones that I can see are in routines that run once per candidate, and they're unlikely to make a big difference.  The code spends two-thirds of its time in function multll(), which multiplies two uint64_t's; I can't see way to improve that routine much.  There may be ways that I can't see, or ways to eliminate calls to it.  

Barrett reduction was something of a chore to implement.  I suspect that it would be a lot easier with BigNumber, but it doing it would require us to get into the raisemod() function in number.c.  I'm not sure I want to try it, and I suspect that it wouldn't beat this speed, anyway.

The code is so long that I can't fit it in with this post.  It's attached, [Edit: attached file, changed this] along with the included file.

Here's output for the representative primes that have become the staple of this thread, from the Uno:
Code:
1000000000054141 is prime. Took 849 mS. Found 4 primes.
1000000000054147 is prime. Took 432 mS. Found 5 primes.
1000000000054259 is prime. Took 995 mS. Found 6 primes.
The Due, of course, ran very fast.  Times were about 10 ms for integers in the high end of the 16-digit range.  Accuracy looks great:  five million primes, from the Due, exactly matched the 46th through 50th million primes from this page: http://primes.utm.edu/lists/small/millions/.

[Edit: Fixed error in attached file.]
124  Community / Exhibition / Gallery / Re: Arduino generating prime numbers on: November 21, 2013, 02:54:25 pm
Thanks.

Wikipedia says so - http://en.wikipedia.org/wiki/Modulo_operation#Equivalencies.  But, in this program, whenever (a*b)%c is calculated, a is either a small prime, or it's <something>%c; and b is either equal to a, or it's <something else>%c.  So, already a%c=a, and b%c=b whenever big numbers are multiplied and modded.  I can't see how to gain ground using that equivalency.

Can you elaborate?
125  Community / Exhibition / Gallery / Re: Arduino generating prime numbers on: November 19, 2013, 07:50:44 pm
One more tweak.  This one doesn't use BigNumber, and it runs nearly twice as fast. 

We used BigNumber to get around the fact that intermediate results can overflow an unsigned long long in the Miller-Rabin test.  That test requires calculating ab % n for large numbers.  Using numbers that are less than 1016 - 54 bits - the results always fit into a long long, but the intermediate products that generate the exponentials can be bigger than 264.  BigNumber gets around the restriction, but at some cost in speed and storage.

This code uses function "mulmodll()," that produces the product of two unsigned long long's, modulo a third unsigned long long, without overflowing on intermediate products.  It generates the product internally as two unsigned long long's, and calculates the modulo by repeated shifts and mods.  It has its own overhead - it has to take <something> % <something else> as many as seven times.  But, it's working in native number formats, rather than in BCD, and it goes faster.

This code depends heavily on the fact that numbers bigger than 1016 aren't investigated.  That fact simplifies the 128-bit multiplication, in that certain intermediate results can't overflow, and it simplifies the mod operations, too.  The code isn't general.  It's assumptions are that the lower half of any argument times the upper half of any argument will be less than 263, and that the modulus will always be less than 256.  With 1016 = 253.15, the problem fits, with a bit or two to spare.

Results:  Running on the Uno, 100000 primes starting at 961,748,929, and comparing to the fiftieth million from this list - http://primes.utm.edu/lists/small/millions/ - found an exact match, with calculation time averaging about 170 ms per prime.  Using the tests of reply #80, calculation time for 16-digit numbers added up to 134,578 ms, for an improvement of a factor of a little better than 5 over using BigNumber, and about 1.9 times as fast as using Numbler4.  For 12-digit numbers, calculation time was 52,472 ms, for a similar speed improvement.  For numbers near the top of the range, starting at 9,999,999,999,990,001, the average calculation time for 100 primes was around 1.5 seconds - about three times what Nick was looking for.

The Due was ferocious.  It averaged 10.6 ms for 1,000 12-digit primes starting at 100000000001; 25.8 ms for 16-digit  starting at 1000000000000001, and 33.4 starting at 9900000000036791.  My Due ran these tests, and then stopped working - I'm hoping that it's just bored.

There's room for improvement here.  I think there's a sweet spot for the number of trial divisions, I think it depends to some extent on the size of the numbers under test, and I don't think I've found it here.  There may be a way to reduce the number of modulo calls in function powmodll(), but I can't think of any yet.

Here's the code.  The included file is attached at reply #70:
Code:
PROGMEM prog_uchar primes[] = {
    3,  5,  7, 11, 13, 17, 19, 23, 29, 31,
   37, 41, 43, 47, 53, 59, 61, 67, 71, 73, 79,
   83, 89, 97,101,103,107,109,113,127,131,137,
  139,149,151,157,163,167,173,179,181,191,193,
  197,199,211,223,227,229,233,239,241,251
};
const byte nPrimes = sizeof(primes)/sizeof(primes[0]);

#include "Pseudoprimes.h";
const uint16_t nPSP = (sizeof(psps)/sizeof(psps[0])) >> 3;
 
uint64_t candidate;

unsigned long start;
long found = 0; // Number of primes found

void setup()
{
  Serial.begin(115200);
  Serial.println ("Starting.");

  candidate = 961748927ULL;
}  // end setup()

void loop()
{
  uint64_t nextPSP;
  int nextPSPindex;
  for (int i=0;i<nPSP;i++) {
    nextPSP = fetchPSP(i);
    if ((nextPSP == candidate) || (nextPSP > candidate)) {
      nextPSPindex = i;
      break;
    }
    nextPSP = 0;
    nextPSPindex = 0;
  }
  start = millis();
 
  for ( ; candidate < 9999999999999999ULL; candidate += 2) {
    if (nextPSP != candidate) {
      if (isPrime(candidate)) {
        found++;
        showPrime (candidate);
      }
    }
    else {
      nextPSPindex++;
      if (nextPSPindex < nPSP) {
        nextPSP = fetchPSP(nextPSPindex);
      }
    }
  }
  found = 1;
  showPrime(1);
  found++;
  showPrime(2);
  candidate = 3;  // back to start
}  // end loop()

bool Miller(uint64_t source)
{
  if(source <= 3) {
    return true;
  }
  if((source & 1) == 0) {
    return false;
  }

  uint64_t d = source - 1;
  int s = 0;

  while((d & 1) == 0)
  {
    d >>= 1;
    s ++;
  }

  uint8_t a;
  const byte bases[] = {2, 3, 5, 7,};// 11, 13, 17, 19, 23 };
  const byte nbases = sizeof(bases)/sizeof(bases[0]);
  for(int i = 0; i < nbases; i++)
  {
    a = bases[i];
    uint64_t x = powmodll(a, d, source);
   if((x == 1) || (x == (source - 1))) {
      continue;
    }

    bool stillOK = false;
    for(int r = 1; r < s; r++)
    {
      x = mulmodll(x, x, source);
      if(x == 1)
        return false;
      if(x == (source - 1)) {
        stillOK = true;
        break;
      }
    }

    if(!stillOK)
      return false;
  }

  return true;
}  // end Miller()

bool isPrime(uint64_t n) {
  if (!check(n)) {
    return false;
  }
  if (n > 0x10000) {
    if (!Miller (n)) {
      return false;
    }
  }
  return true;
}

void showPrime (uint64_t n)
{
  printll(n);
  Serial.print(" is prime. Took ");
  Serial.print(millis () - start);
  Serial.print (" mS. Found ");
  Serial.print (found);
  Serial.println (" primes.");
  start = millis ();

}  // end showIfPrime()

bool check(uint64_t n) {
  if (n <= 3) {
    return true;
  }
  if (!(n & 1)) {
    return false;
  }
  uint8_t i = 1;
  uint8_t f = pgm_read_byte_near(primes);
  while ((i < nPrimes) && (f < n)) {
    if ((n % f) == 0) {
      return false;
    }
    f = pgm_read_byte_near(primes+i);
    i++;
  }
  return true;
} // end check()

uint64_t fetchPSP(int i) {
  i <<= 3;
  uint64_t n = 0;
  for (int8_t k=0;k<sizeof(n);k++) {
    n = n <<= 8;
    n |= pgm_read_byte_near(psps+i+k);
  }
  return n;
} // end fetchPSP()

void printll(uint64_t n) {
  // print an unsigned long long
  if (n == 0) {
    Serial.print('0');
    return;
  }
  const byte DIGITS = 21;
  char s[DIGITS+1];
  for (byte i=0;i<DIGITS;i++) {
    s[i] = 'X';
  }
  s[DIGITS] = 0;
  byte i = DIGITS-1;
  while (n > 0) {
    s[i] = (n % 10) | '0';
    n = n / 10;
    i--;
  }
  i++;
  char *z = &s[i];
  Serial.print(z);
} // end printll()

uint64_t powmodll(uint8_t a, uint64_t b, uint64_t c) {
  // Calculate a**b % c,
  uint64_t x = 1;
  uint64_t y = a;

  while(b > 0) {
    if(b & 1) {
       x = mulmodll(x, y, c);
      }
      y = mulmodll(y, y, c);
      b >>= 1;
    }
    return x;
} // end powmodll()

uint64_t mulmodll(uint64_t a, uint64_t b, uint64_t c) {
  // Calculate a*b %c
  uint64_t imh; // Intermediate results for multiplication
  uint64_t iml;
  uint64_t imm;

  // Multiply a and b
  // Routine is not general:
  //   relies on fact that MSB of inputs is zero
  // Product stored as two unsigned long long's
  imm = (a & 0xffffffffUL) * (b >> 32);
  iml = (a >> 32) * (b & 0xffffffffUL);
  imm += iml; // overflow for general inputs
  iml = (a & 0xffffffffUL) * (b & 0xffffffffUL);
  imm += iml >> 32;
  iml = (iml & 0xffffffffUL) + (imm << 32);
  imh = (a >> 32) * (b >> 32);
  imh += imm >> 32;

  if (imh == 0) {
    return (iml % c);
  }
  for (uint8_t i=0;i<8;i++) {
    imh = imh << 8;
    imh |= iml >> 56;
    iml = iml << 8;
    if ((imh >> 56) != 0) {
      if (imh > c) {
        imh %= c;
      }
      else {
        if (imh == c) {
          imh = 0;
        }
      }
    }
  }
  if (imh > c) {
    imh %= c;
  }
  return imh;
} // end mulmodll()

Here's Uno output for comparison to reply #72, with clock speed the standard 16 MHz:
Code:
1000000000054141 is prime. Took 1433 mS. Found 1 primes.
1000000000054147 is prime. Took 778 mS. Found 2 primes.
1000000000054259 is prime. Took 1594 mS. Found 3 primes.
126  Community / Exhibition / Gallery / Re: Arduino generating prime numbers on: November 13, 2013, 01:41:41 am
I don't understand why you seem to think that BigNumber is so quick.
No offense intended.  BigNumber is quick compared to my expectations.  I'm surprised to see how well it does compared to operating with native integers.  I had expected it to have significant overhead.  

OK, here are results with BigNumber, and with Numbler4 as posted at http://forum.arduino.cc/index.php?topic=85692.msg1416890#msg1416890.  I transplanted powMod() from BigNumber into Numbler4.  Tests ran on an Uno and a Due, and, just for fun, with a TI Stellaris Launchpad, compiled under the open-source Arduino-like IDE, "Energia. " Calculated primes agreed for all three devices.  Results are total time in milliseconds for calculation of the first 100 primes after 10N-1, where N is 12 and 16:

                         12 digit                16 digit
           Uno:      102002    263207        254004    688223
           Due:       14431     13790         34571     32051
     Stellaris:       10349     10281         24643     24058


It ran significantly faster with Numbler4 on the Uno, with execution time less than 40% of the BigNumber version.  It ran a little slower than BigNumber on the ARMs.  On the Due, Numbler4 ran about 5-8% slower.  There was less difference with the Stellaris, but it was still bit slower.  I can qualitatively understand how eliminating a lot of modulo operations can speed things up on the Uno; my understanding of the internals of the ARMs is too thin to help me explain the results with the Due and Stellaris.  I'm fool enough to speculate, though:  TI says that the Stellaris' Cortex-M3 has hardware divide, and Atmel says the same thing about the SAM3X.  My guess is that division on those processors is fast enough to be competitive with the alternative.

Uno and Due programs were identical, generally the code from reply 66, maybe with a minor change or two.  The unused function rng(), for generating random bases for the Miller-Rabin test, was deleted.  The IDE v1.5.4 wouldn't compile it for the Due - it didn't like the random() function without arguments, preferring rand() instead.  v1.5.2 was OK with random().  I don't even want to think about about the machinations I had to go through to get the Stellaris to run the program; just as well, since that kind of talk doesn't belong on an Arduino forum.
127  Community / Exhibition / Gallery / Re: Arduino generating prime numbers on: November 11, 2013, 11:26:17 pm
One last tweak.  This one includes a list of numbers that fool the Miller-Rabin test for bases 2, 3, 5, and 7.  It compares each odd number to the next item in the list, and skips it if it finds them equal.  The previous version, in reply #66, uses nine bases, as required to eliminate pseudoprimes below 1016. This one uses four and a table lookup, so it's a lot faster.  It kind of feels like cheating, though - Nick's versions knew that 2 an 3 are prime, while this program has a lot of arcane knowledge about special numbers.  

I was surprised to see BigNumber work so quickly compared to native number formats, as I'd expect it to carry substantial overhead.  I suspect that its progenitors were developed, at least in part, for dealing with calculations of primes, and that may account for its fitness for this use.

My thanks to Nick for posting this project, and for coming up with the BigNumber library. Hunting for primes is easy; speeding up the process is really tricky.  And, I learned how to spell, "pseudo-," I think, so my ability to eruditely describe things that aren't quite what they seem is much improved.

A list of strong psuedo fake primes to bases 2 and 3 is here - http://www.bell-labs.com/user/bleichen/diss/thesis.html.  The CPAN module Math::Primality provided the engine for extracting the pseudoprimes that penetrate tests for bases 5 and 7.  The included file containing the byte-wide bits of the list is attached - "Pseudoprimes.h".

The code has its oddities.  Because there seems to be no way to store or access uint64_t values with PROGMEM, and because I couldn't get BigNumber to play nice with variables bigger than one byte, pseudoprimes are stored as 8 bytes, big-endian, and decoded into BigNumbers.  Here's the sketch:
Code:
#include <BigNumber.h>
PROGMEM prog_uchar primes[] = {
    2,  3,  5,  7, 11, 13, 17, 19, 23, 29, 31,
   37, 41, 43, 47, 53, 59, 61, 67, 71, 73, 79,
   83, 89, 97,101,103,107,109,113,127,131,137,
  139,149,151,157,163,167,173,179,181,191,193,
  197,199,211,223,227,229,233,239,241,251
};
const byte nPrimes = sizeof(primes)/sizeof(primes[0]);

#include "Pseudoprimes.h";
const uint16_t nPSP = (sizeof(psps)/sizeof(psps[0])) >> 3;
  
BigNumber candidate;
BigNumber one (1);
BigNumber two (2);

unsigned long start;
long found = 0; // Number of primes found

void setup()
{
  Serial.begin(115200);
  Serial.println ("Starting.");

  BigNumber::begin ();  

  candidate = BigNumber ("100000000003");
}  // end of setup

void loop()
{
  BigNumber nextPSP;
  int nextPSPindex;
  for (int i=0;i<nPSP;i++) {
    nextPSP = fetchPSP(i);
    if ((nextPSP == candidate) || (nextPSP > candidate)) {
      nextPSPindex = i;
      break;
    }
    nextPSP = 0;
    nextPSPindex = 0;
  }

  start = millis();
  
  for ( ; candidate < BigNumber("9999999999999999"); candidate += 2) {
    if (nextPSP != candidate) {
      showIfPrime (candidate);
    }
    else {
      nextPSPindex++;
      if (nextPSPindex < nPSP) {
        nextPSP = fetchPSP(nextPSPindex);
      }
    }
  }
  candidate = 3;  // back to start!
}  // end of loop

bool Miller(BigNumber source)
{

  if(source == 2 || source == 3)
    return true;
  if(source < two || source % two == 0)
    return false;

  BigNumber d = source - one;
  int s = 0;

  while(d % two == 0)
  {
    d /= two;
    s += one;
  }

  BigNumber a;
  const byte bases[] = {2, 3, 5, 7,};
  const byte nbases = sizeof(bases)/sizeof(bases[0]);
  for(int i = 0; i < nbases; i++)
  {
    a = BigNumber (bases[i]);
    BigNumber x = a.powMod (d, source);
    if(x == one || x == source - one)
      continue;

    bool stillOK = false;
    for(int r = 1; r < s; r++)
    {
      x = x.powMod(two, source);
      if(x == one)
        return false;
      if(x == source - one) {
        stillOK = true;
        break;
      }
    }

    if(!stillOK)
      return false;
  }

  return true;
}  // end of Miller

void showIfPrime (BigNumber n)
{
  if (!check(n))
    return;
  if (!Miller (n))
    return;
  found++;
  Serial.print(n);
  Serial.print(" is prime. Took ");
  Serial.print(millis () - start);
  Serial.print (" mS. Found ");
  Serial.print (found);
  Serial.println (" primes.");
  start = millis ();

}  // end of showIfPrime

bool check(BigNumber n) {
  BigNumber f;
  for (byte i = 0;i < nPrimes; i++) {
    f = BigNumber(pgm_read_byte_near(primes+i));
    if ((n % f) == 0) {
      return false;
    }
  }
  return true;
}

BigNumber fetchPSP(int i) {
  // Progmem array psps contains strong pseudoprimes
  // 8-byte numbers, arranged as bytes, big-endian
  i <<= 3;
  BigNumber x256 ("256");
  BigNumber n = 0;
  for (byte k=0;k<8;k++) {
    n = n * x256;
    uint8_t b = pgm_read_byte_near(psps+i+k);
    n = n + BigNumber(b);
  }
  return n;
}
Here's some output:
Code:
100000000003 is prime. Took 1676 mS. Found 1 primes.
100000000019 is prime. Took 2798 mS. Found 2 primes.
100000000057 is prime. Took 3539 mS. Found 3 primes.
100000000063 is prime. Took 1770 mS. Found 4 primes.
Here's the output for a Due running the same code:
Code:
100000000003 is prime. Took 86 mS. Found 1 primes.
100000000019 is prime. Took 145 mS. Found 2 primes.
100000000057 is prime. Took 189 mS. Found 3 primes.
100000000063 is prime. Took 89 mS. Found 4 primes.
[Edit: Add this] And, the Due in the 16-digit zone:
Code:
1000000000000037 is prime. Took 416 mS. Found 1 primes.
1000000000000091 is prime. Took 236 mS. Found 2 primes.
1000000000000159 is prime. Took 529 mS. Found 3 primes.
1000000000000187 is prime. Took 261 mS. Found 4 primes.
128  Using Arduino / Audio / Re: Amplified Signal to Arduino Mega on: November 07, 2013, 05:22:42 pm
... automatic gain control amplifiers, I assume this IC's advertised automatic level control is the same thing?
Same-ish.  Automatic level control, as used in the datasheet, implies that the output will be within a vary narrow range.  

Quote
"...voltage input range of 1000." Do you mean over a range of 1000dB, 1000V, 1000x, or what?
A range of 1000x.  If the datasheet's claims are accurate, you should be able to design a circuit with that will deliver a constant AC voltage at the output while the input signal varies by a factor of 1000.

Quote
I have drawn a schematic based on the diagram in the data sheet and attached. Does it look correct?
I don't know.  I'm not an expert in using the SA575.  I have an application for this IC myself, and I was hoping that you'd design and debug something, and then tell me how to use it.  I think that you're in for some experimentation.  My suggestion is to implement the circuit shown in the datasheet, and see how it performs.  Armed with that information, you'll be able to calculate how to manage the IC's output to use it as input to something else.

I'd recommend considering the output voltage of your amp to be more like sqrt(2) * sqrt(50*6), and then add a safety factor, and then protect your inputs to avoid damage from overvoltage.  The output power quoted for the amplifier is probably imprecise, and it's based on calculations that might be different from the ones that you and I would make.  After you know more, then you can hook it up to the next delicate circuit.

129  Community / Exhibition / Gallery / Re: Arduino generating prime numbers on: November 05, 2013, 01:14:09 pm
Got a speed improvement of an average factor of about 2.3 over the Nick's BigNumber-based code shown in reply #43, starting at 1000000000000001, 16 digits, by making these changes:
  • Added a trial division test using byte-sized primes, before other testing.
  • Modified the Miller-Rabin test to use 2, 3, 5, 7, 11, 13, 17, 19, and 23, rather than random numbers.  This Wikipedia link - http://en.wikipedia.org/wiki/Strong_pseudoprime#Examples - says that there's no number less than 3825123056546413051, 19 digits, about 1.6 x 261, that will fool every single-pass Miller-Rabin tests with these bases.
  • [Edit: Add this] Modified loop parameters to wield only integers = 6k + 1 and 6k - 1, for integer k.

Most of the speed improvement appears to come from eliminating Miller-Rabin tests for numbers that divide by small primes.  That suggests that there may be a sweet spot for the number of prime divisors in the trial division test that depends on the speed of the Miller-Rabin test, and it may be different from the one I chose.  Lesser improvement comes from eliminating selection of random BigNumbers, and reducing the number of Miller-Rabin tests from ten to nine.  [Edit: Add this] The 6K +- 1 trick doesn't do much, since it only eliminates multiples of 3, and those fail the trial division test almost immediately.

Average displayed time was about 11.5 mS for the primes shown in replies #44 and #45.   It also ran about 3.5% slower than the values shown in reply #64 - calculation times were sometimes wildly different, but on average, they were almost as good.

Here's the code, based on the sketch in reply #43, minimally changed, auto-formatted:
Code:
#include <BigNumber.h>
PROGMEM prog_uchar primes[] = {
  2,  3,  5,  7, 11, 13, 17, 19, 23, 29, 31,
  37, 41, 43, 47, 53, 59, 61, 67, 71, 73, 79,
  83, 89, 97,101,103,107,109,113,127,131,137,
  139,149,151,157,163,167,173,179,181,191,193,
  197,199,211,223,227,229,233,239,241,251
};
const byte nPrimes = sizeof(primes)/sizeof(primes[0]);

const char FIRST_PRIME_TO_USE [] =  "1000000000000001";

BigNumber candidate;
BigNumber one (1);
BigNumber two (2);

byte bump = 2;

void setup()
{
  BigNumber five (5);
  BigNumber six (6);

  Serial.begin(115200);
  while (!Serial) {
  }

  BigNumber::begin ();  

  Serial.println ();
  Serial.println ();
  Serial.println ("Starting.");

  candidate = BigNumber (FIRST_PRIME_TO_USE);

  BigNumber rem = candidate % six;
  while ((rem != one) && (rem != five)) {
    rem = rem + one;
    rem = rem % six;
    candidate = candidate + one;
  }
  if (rem == one) {
    bump = 4;
  }
}  // end of setup

void loop()
{

  for ( ; candidate < BigNumber ("9999999999999999"); candidate += bump, bump = 6 - bump)
    showIfPrime (candidate);

  candidate = 3;  // back to start!
}  // end of loop

void rng (BigNumber & result)
{
  result = BigNumber (random ());  
}  // end of rng

bool Miller(BigNumber source)
{

  if(source == 2 || source == 3)
    return true;
  if(source < two || source % two == 0)
    return false;

  BigNumber d = source - one;
  int s = 0;

  while(d % two == 0)
  {
    d /= two;
    s += one;
  }

  const byte bases[] = { 2, 3, 5, 7, 11, 13, 17, 19, 23 };
  const byte nbases = sizeof(bases)/sizeof(bases[0]);
  for(int i = 0; i < nbases; i++)
  {
    //    do
    //    {
    //      rng (a);
    //    }
    //    while(a < two || a >= source - two);

    BigNumber a = BigNumber (bases[i]);
    BigNumber x = a.powMod (d, source);
    if(x == one || x == source - one)
      continue;

    for(int r = 1; r < s; r++)
    {
      x = x.powMod(two, source);
      if(x == one)
        return false;
      if(x == source - one)
        break;
    }

    if(x != source - one)
      return false;
  }

  return true;
}  // end of Miller

const byte DIGITS = 16;

unsigned long start;
unsigned long lastShown;
long found = 1; // Number we found

void showIfPrime (BigNumber num)
{
  if (!check(num))
    return;
  if (!Miller (num))
    return;

  char buf [20] = { 0 };
  long long temp = num;
  byte pos = 0;

  // make printable
  for (pos = 0; pos < DIGITS; pos++)
  {
    if (temp == 0)
      break;
    char digit = temp % 10;
    buf [DIGITS - 1 - pos] = digit | '0';
    temp /= 10;
  }

  // insert leading spaces
  for (; pos < DIGITS; pos++)
    buf [DIGITS - 1 - pos] = ' '; //  = 15;


  Serial.print(num);
  Serial.print(" is prime. Took ");
  Serial.print(millis () - start);
  Serial.print (" mS. Found ");
  Serial.print (found);
  Serial.println (" primes.");
  start = millis ();

  found++;

  lastShown = millis ();

}  // end of showIfPrime

bool check(BigNumber n) {
  for (byte i = 0;i < nPrimes; i++) {
    BigNumber f = BigNumber(pgm_read_byte_near(primes + i));
    if ((n % f) == 0) {
      return false;
    }
  }
  return true;
}
130  Using Arduino / Audio / Re: Amplified Signal to Arduino Mega on: October 31, 2013, 09:01:37 pm
This gizmo claims that it can operate as a limiter, and accommodate 60 dB on its inputs with only about 0.5 dB change on the outputs. 

http://www.onsemi.com//pub/Collateral/SA575-D.PDF

It's available in a through-hole mount for something around $5.50 from Mouser, though the through-hole version is marked "end-of-life."  If that's true, you should be able to get the a flat output response over a voltage input range of 1000.  Using the "+10 dB ~ twice as loud" rule of thumb, that's an apparent loudness range factor of 64.  That might do it.
131  Using Arduino / Programming Questions / Re: kurtosis calculations on an analog signal on: October 29, 2013, 07:32:50 pm
Every time loop() executes, the counter variable - i - is incremented, whether it's time to take a sample or not.  loop() executes many times between samples, so the counter reaches its limit quite quickly, and the output data prints.

Since the counter is supposed to be counting samples, it needs to be incremented only when a sample is taken.  It needs to be inside the brackets of the "if (currentMillis ... " test.

[Edit: Add this]  After a closer look:
  • Proper indentation will make your code easier to read and understand, for others and for you.  You can manage that yourself, or you can use "Auto Format" under the "Tools" menu in the Arduino IDE.
  • I recommend that you refrain from adding complexity to this code until you get the kurtosis calculation working.  Leave xbee stuff out of your program until you get kurtosis working.
  • Patience, Grasshopper.
132  Using Arduino / Programming Questions / Re: kurtosis calculations on an analog signal on: October 29, 2013, 11:39:58 am
Don't forget:  the first goal is to get the calculation working.  It's a whole lot easier to modify a working program than to start over with something new.  Fix this one, and worry about optimization and and improvements later.
133  Using Arduino / Programming Questions / Re: kurtosis calculations on an analog signal on: October 29, 2013, 08:15:10 am
Ingenious.  That's a lot of reasonable-looking code to come from a guy who didn't know where to start a few days ago.

Here are some suggestions:
  • [Edit:  Add this]
    Something's wrong with the timing and program flow.  The code decides when to take a sample, but it processes a sample regardless of whether or not it takes one.  Examine it:  loops starts; the test passes, it takes a sample, and processes it; loop terminates, and starts over; the test fails, and it doesn't take a sample, but it processes the last sample again anyway.  Most of that code needs to move into the "if(currentMillis ..." section, and the rest should only be executed when it's time to do final calculations and print.
  • You're trying to calculate a function of a series of values, and you don't yet know whether the calculation is correct.  You also don't know what the values are, since they're analog readings, and they aren't displayed.  Troubleshooting the calculation is difficult enough without having to wonder what the input data was.  While you're troubleshooting the calculation, use some convenient calculated data whose kurtosis is well-known, instead of analog readings.  Shutting off sample timing during troubleshooting will make the results show up faster.  Here are two suggestions:
Code:
xi = i  % 10
          // uniform distribution, expected excess kurtosis = -1.2
xi = sin(i*0.001*2*3.1415927*100) + 1
         // sine function, 100 Hz, sampled at 1/millisecond, expected excess kurtosis = -1.5

  • Decide whether you want to calculate kurtosis the old-fashioned way, or whether you want to calculate excess kurtosis.
  • Print some debugging information, like the value of the running average, fourth moment, and variance.  Be prepared to get deeper into the calculations if those values don't help.
  • Your "running average" technique, where newValue = alpha*sampleValue + (1-alpha)*oldValue, is a single-pole digital low-pass filter.  Its output value tends toward the mean of the input signal.  When you filter variance, you get something that looks like the statistical variance - the mean of the squares of the differences between the samples and the mean of the signal.  When you filter the fourth moment, you get something that looks like the fourth moment - the mean of the fourth power of those same differences.  Division by the number of samples is essentially accomplished in filtering those values.   Consequently, I'll recommend that you reconsider whether this line needs the multiplication by the number of samples:
Code:
 denominator = sigma4*numberSamples;
  • When you start getting output data, experiment with larger values of alpha.  I don't think it's necessary that alpha be the same as 1.0/numberSamples.  You might want to define beta  = 1.0 - alpha.
  • This looks to me to be an ingenious solution.  My experience with ingenious solutions to statistical calculations is that they often yield unpextected results, usually because there's some nugget of probability theory that I've overlooked that makes them invalid.  I can't see the flaw here, though - the level of filtering is pretty intense, and the signal mean will get pretty close to the actual value; once it stabilizes, I think that the variance and fourth moment calculations will be valid.  Nevertheless, I'd recommend testing it with several kinds of distributions to see that it gives the expected results, and be prepared for disappointment.
  • After you're confident of the kurtosis calculation, think about ways to mange the calculation of the signal mean.  When I run this code with known input data and print the value of runningaverage, it takes a while to stabilize.  Until it does, the output data will not be valid.  If the DC offset of the signal won't change during the measurement, consider instead taking some time to determine the signal mean while outputting no values, and then using that mean to calculate the second and fourth moments.  There's no loss of performance - there's a delay between startup and valid data either way.  With the mean known, the calculations are pretty simple.  Once the signal mean is established, I think that there's no necessity to filter the second and fourth moments - the calculation can be straightforward from the formula.
  • Consider AC-coupling the input signal.  The signal mean will then be the value of the analog offset voltage.  That'll give control of the signal mean to you, rather than to some setting on a sensor amplifier.  There are a lot of options for how to manage getting that voltage into the calculation, some involving hardware, some not.  A simple one would be to use a resistor divider to set it at half scale, and use that value as the initial value of runningaverage, rather than starting with zero.  That would reduce the time it takes for data to become valid, though possibly not by all that much - the second and fourth moments are still filtered.
  • As PaulS suggests, don't make sensorPin a float.  It seems to work, but it's just too weird.  sensorPin doesn't change, so it could easily be const int.
134  Using Arduino / Programming Questions / Re: kurtosis calculations on an analog signal on: October 28, 2013, 04:25:12 pm
It's good to see the original poster put up some code.  It suggests that he's actually interested in a solution.

The posted code doesn't compile.  There are repeated references to kurtvoltages, as if it were a variable. But kurtvoltages is an array, and, if it's going to be treated like a single value, it needs an index.

There are logic errors, too, but they can wait.  If it doesn't compile, you'll never get any results at all.

[Edit: Add this]
And, until you know what you're doing, please don't try to hook up an analog source - like a signal generator or a biological monitor - to your Arduino.  It can be quite easy to damage an analog input when the input voltage goes out of bounds.  I suspect that your signal source is centered at zero volts, and that it's negative for at least some of the time.  A negative voltage on an analog pin can damage it.

This sounds like an academic assignment.  Is it?
135  Using Arduino / Project Guidance / Re: Measuring AC line frequency on: October 08, 2013, 05:10:34 am
Here are three application notes from Atmel that might help:


Pages: 1 ... 7 8 [9] 10 11 ... 25