[alex1985 0]

Здравствуйте, уважаемые форумчане!

 

Задача в следующем: есть таймер, например, 5-й, питается от 72 мГц шины. Требуется генерировать ШИМ частотой в пределах от 10 до 2000 Гц с минимально возможной погрешностью (разницы нет - в большую или меньшую сторону). Когда частота кратна 72 - все прекрасно, но если требуется вывести частоту, например, 13 Гц, то возникает погрешность. Настройка соответственно происходит регистрами PSC и ARR, комбинаций их может быть достаточно много. Самое первое, что приходит в голову - это тупой перебор всех возможных комбинаций с вычислением ошибки и выбором пары PSC и ARR с наименьшей их них, но это не наш метод ;)

Может существует формула для расчета? Или на крайний случай таблица значений регистров для каждой частоты с минимальной погрешностью?

[Eddy_Em 1]

Подсказываю: есть цепные дроби. Я, кстати, активно их использую для получения коэффициентов преобразования ADU с АЦП в единицы измерения. В octave есть функция rat, которая эти самые цепные дроби считает. Скажем, надо аппроксимировать число 0.13, пишем rat(0.13) и получаем:

octave:2> rat (0.13)
ans =

                           
0 + 1/(8 + 1/(-3 + 1/(-4)))

Итак, простейшее приближение - это 1/8 (3.8% погрешности). Следующее — 1/(8-1/3)=1/(23/3)=3/23 (погрешность 0.3%). Последнее приближение даст нам уже 1/(8 - 4/13)=1/(100/13)=13/100 (логично, да? ☺)

И так с любым дробным числом.

 

Выбираем дробь, удовлетворяющую заявленной погрешности и соответствующим образом настраиваем PSC и ARR.

 

Кстати, у STM32 есть еще и возможность связывания таймеров. Так что, можно очень точно генерировать такие низкие частоты. Вот в мегагерцевом диапазоне - да, проблемка.

[adnega 10]

есть таймер, например, 5-й, питается от 72 мГц шины. Требуется генерировать ШИМ частотой в пределах от 10 до 2000 Гц с минимально возможной погрешностью

Вопрос понятен, но не легче ли использовать PSC=0 и 32-битный ARR ? TIM2 в помощь ))

[BaN 0]

Подсказываю: есть цепные дроби.

Смотря какую "минимально возможную погрешность" ищет ТС, я так полагаю, что не требуется точность до сотых долей Гц, так что если взять частный случай с пределами от 10 до 2000 Гц, то дроби можно откидывать сразу.

Fout = Fclk / (PSC +1) / (ARR + 1)

(PSC +1) * (ARR + 1) = Fclk / Fout

Если взять худший случай, что Fout = 2000 Гц, то получим (Fclk / Fout) = 72000000 / 2000 = 36000 и здесь уже нет смысла учитывать дробную часть.

Задача сводится к разложению целого числа (Fclk / Fout) на целые множители с учетом пределов изменения множителей 1...65536.

[_Pasha 0]

В общем-то решается не только цепными дробями.

Я тут тоже посидел, по случаю, накрапал вот:

#include <stdint.h>
#include <stdlib.h>
#include <stdio.h>

static unsigned int fast_sqrt(unsigned int arg)
{
   unsigned int one,res = 0;
   for(one = 1UL << (sizeof(one)*8 - 2); one > arg; one >>=2);

   while(one)
   {
      unsigned int k= res+one;
      if(arg >= k)
      {
           arg -= k;
           res += 2*one;
      }
      res >>= 1;
      one >>= 2;
   }
   return res;
}

unsigned int factorize(unsigned int arg, unsigned int initial)
{
   unsigned int j, limit = fast_sqrt(arg);
   if(initial == 1) initial += 1;
   for(j=initial; (j<=limit)&&(arg%j);++j);
   return (j<=limit)?j:0;
}

void get_mpys(uint16_t *a, uint16_t *b, unsigned int num)
{
unsigned int koeff=1;
*a = *b = 1;
while(koeff)
{
	unsigned int k = factorize(num,koeff);
	if(k)
	{
		koeff = k;
		num /= k;
	}
	else
	{
		k = num;
		koeff = 0;
	}

	if(*a * k < 65536)
		*a *= k;
	else
	{
		if(*b * k < 65536)
			*b *= k;
		else
		{
			*a=*b=0;
			return;
		}
	}
}
}

/*test*/
int main()
{
uint16_t xa,xb;
unsigned int nnn= 7200000/7;
get_mpys(&xa,&xb, nnn);
printf("numb %d  a= %d  b= %d  a*b = %d\n",nnn,xa,xb,xa*xb);
return 0;
}

 

Вкратце. Имеем static unsigned int fast_sqrt(unsigned int arg) - быстрый целочисленный квадратный корень . Он нам нужен для оценки верхнего предела факторизации, то бишь разложения на множители.

Далее, сам факторизатор - unsigned int factorize(unsigned int arg, unsigned int initial)

Работает тупым перебором, находит множитель начиная с числа initial для параметра arg

Наконец, void get_mpys(uint16_t *a, uint16_t *b, unsigned int num) - тут уже идет последовательное нахождение коэффициентов и домножение их к a или b

То есть a*b=num

Если число не помещается в 16 разрядов - пардон, пользуемся ARR 32 бита. :smile3046:

Можно сравнить погрешности если кто цепные дроби подымет здесь для тестинга.

Изменено 17 мая, 2015 пользователем _Pasha

[basileus_ 0]

предложу следующее решение:

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

const unsigned int FQz=24000000;//входная частота таймера.
const float erratum=.00001;//допустимая ошибка.
const unsigned int MAXLOAD=65535;
//unsigned int  Fs; //требуемая для установки частота.


void calc(unsigned int Fs, unsigned int *ld1,  unsigned int *ld2){
int delta1;
int delta2;
unsigned int ldt1;
unsigned int ldt2;
unsigned int mult=(unsigned int)floor(FQz/(Fs) + 0.5);// получим общее произведение делителя частоты
unsigned int ldx=mult/65535; if(ldx==0){++ldx;}; //возьмём нижнюю оценку для делителя 1.

*ld2=0;
delta2=FQz;// загрузим заведомо превыщающее оценку значение
int cnt=0;
while((ldx<sqrt(mult)))
{
ldt2=(unsigned int)floor( mult/ldx);
ldt1=(unsigned int)floor(mult/ldt2);
if ((ldt2<MAXLOAD)){
	delta1=abs(FQz-Fs*ldt1*ldt2);
	if(delta1>abs(FQz-Fs*ldt1*(ldt2+1))){
	++ldt2;
	delta1=abs(FQz-Fs*ldt1*ldt2);
	}
	if ((delta1<delta2)){//обновим пару делителей, при вновь найденной минимальной ошибке. 
		*ld1=ldt1;
		*ld2=ldt2;
		delta2=delta1;
		};
	if(delta1==0)break;//если найден идеальный результат- остановимся.
	//if((float)delta1<erratum*FQz)break;//если найден подходящий результат- остановимся.
	if((delta1<erratum*FQz)&&(cnt==0)){cnt=ldx;}//если найден подходящий результат- запомним.
};
++ldx;
};
printf("\nfs=%d;  ld1=%u; ld2=%u; delta=%d; ldx=%i; cnt=%i", Fs, *ld1, *ld2, delta2, ldx, cnt );//чисто отладочный вывод
}

int main() {
unsigned int ld1;
unsigned int ld2;

for (int i=1; i<999;i++){
calc(i, &ld1, &ld2);
};
printf("This is SPARTA!\r\n");
}

можно убрать убрать строки, где считается число итераций (cnt), и сразу по достижении приемлемого результата брякаться.

 

 

[adnega 10]

Как вариант, для низких частот (когда делитель не помещается в 16 бит) генерировать циклы с разными ARR.

Например, для 11Гц делитель получается 6545454 = 99 * 65536 + 57390. Устанавливаем PSC = 0, ARR = 65536 - 1, и в обработчике прерывания по UF

считаем до 99 после этого устанавливаем ARR = 57390 - 1. При возникновении прерывания по UF устанавливаем ARR = 65536 - 1. Все повторяется.

Единственный момент, нужно следить, чтоб остаток не приближался к 65535. В этом случае можно в конце сделать два цикла 1/4 остатка и 3/4 остатка:

14347 и 43043 соответственно. Гарантируется, что при возникновении прерываний будет достаточно времени установить ARR до момента прохода таймером

этого значения. Кста, можно в разных циклах накапливать ошибку и получить точные 11Гц за счет того, что иногда цикл на 1 такт длиннее.

 

UPD: Нагрузка на CPU чуть-чуть повышается, но ошибку генерации частоты можно получить нулевую.

[alex1985 0]

Спасибо всем за советы! Сейчас буду перемалывать )

[jcxz 184]

Я тут тоже посидел, по случаю, накрапал вот:

...

unsigned int factorize(unsigned int arg, unsigned int initial)

{

unsigned int j, limit = fast_sqrt(arg);

if(initial == 1) initial += 1;

for(j=initial; (j<=limit)&&(arg%j);++j);

return (j<=limit)?j:0;

}

Это у Вас разложение на простые множители?

При разложении большого arg будет много лишних операций деления.

Возможно лучше после оценки предела (посредством sqrt()), найти все простые числа до него с помощью "решета Сундарама" и пытаться делить только на них.

 

PS: Кстати - в какую сторону даёт погрешность fast_sqrt()? Если в меньшую - надо быть осторожным при разложении arg == n * n, где n - простое число.

Увеличить бы предел на неск. процентов: limit = fast_sqrt(arg) * 11 / 10 или использовать: limit = 1 << (32 - __CLZ(arg) + 1) / 2

[_Pasha 0]

PS: Кстати - в какую сторону даёт погрешность fast_sqrt()? Если в меньшую - надо быть осторожным при разложении arg == n * n, где n - простое число.

кстати да. нужна финальная проверка.

 

 

По первому вопросу не вникал. Яснопонятно, что это все дороговатые операции, однако, таких приложений сферических , чтобы генерить произвольную частоту - еще поискать нужно.

Проще таблицу простых чисел заиметь во флеше, например

http://integraloff.net/simple/table_simpl_ch.php

до 10**4 это всего-то 2460 байт

 

unsigned int factorize(unsigned int arg, unsigned int initial)
{
unsigned int j, limit = fast_sqrt(arg);
if(initial == 1) initial += 1;
if(limit*limit < arg) limit +=1;
for(j=initial; (j<=limit)&&(arg%j);++j);
return (j<=limit)?j:0;
}

[jcxz 184]

Проще таблицу простых чисел заиметь во флеше, например

http://integraloff.net/simple/table_simpl_ch.php

до 10**4 это всего-то 2460 байт

Для решета Сундарама нужно всего 4килобайта ОЗУ для поиска всех простых чисел в диапазоне до 65536. А этого достаточно для разложения на множители любого 32-битного числа.

И алгоритм использует только умножения-сложения.

[_Pasha 0]

Для решета Сундарама нужно всего 4килобайта ОЗУ для поиска всех простых чисел в диапазоне до 65536. А этого достаточно для разложения на множители любого 32-битного числа.

И алгоритм использует только умножения-сложения.

Интересно. А у нас же ж факторизуемое число уменьшается и нижний предел поиска натуральных чисел растет.

Получается, что и столько озу ему может быть не нужно.

[adnega 10]

А этого достаточно для разложения на множители любого 32-битного числа.

Если это число не будет простым.

[_Pasha 0]

Для решета Сундарама нужно всего 4килобайта ОЗУ для поиска всех простых чисел в диапазоне до 65536.

:yeah: Вчера мы не сообразили, что есть же еще компромисс, а именно:

если таблицу всех простых чисел представить битовым полем, то до 65536 нужно 8к флеша. и это без RLE.

А функция распределения простых числел https://ru.wikipedia.org/wiki/%D0%A4%D1%83%...%81%D0%B5%D0%BB

так и шепчет его применить.

 

upd

 

def prime_nums(N):
    return reduce( (lambda r,x: (r.difference_update(range(x*x,N,2*x)) or r)
        if (x in r) else r), range(3, int((N+1)**0.5+1), 2),
        set([2] + range(3,N,2)))
        
aList= prime_nums(256**2)
nprims=len(aList)
print "//found {0} prime numbers".format(nprims)
#create u32-bitmap
u32bmp=[]
u32bmp[:]= [0 for j in xrange(65536/32)]
#map to
for nn in aList:
    u32bmp[nn/32] |= 1 << (nn % 32)
print "static const uint32_t prime_set[]={"
col = 0
ss=''
for nn in u32bmp:
    print "0x%08X," % nn
print "};"

8кб

дальше я не соображаю, как его топтать.

Может кто напишет?

Или хотя бы проверит, правильно ли работает?

upd2 -поменял вывод, а то некрасиво

Изменено 22 мая, 2015 пользователем _Pasha

[_Pasha 0]

вишенка на торте

char isPrime(uint_least16_t num)
{
    return (prime_set[num/32] & (1<<(num & 31)))? 1:0;
}

unsigned int next_prim(unsigned int current)
{
   while(!isPrime(++current));
   return current;
}

unsigned int factorize(unsigned int arg, unsigned int initial)
{
  unsigned int j, limit = fast_sqrt(arg);
  if(initial == 1) initial += 1;
  if(limit*limit < arg) limit +=1;
  for(j=initial; (j<=limit)&&(arg%j); j=next_prim(j));
  return (j<=limit)?j:0;
}