diff options
Diffstat (limited to 'mat/programčki/ničle.singlethreaded.c')
-rw-r--r-- | mat/programčki/ničle.singlethreaded.c | 74 |
1 files changed, 74 insertions, 0 deletions
diff --git a/mat/programčki/ničle.singlethreaded.c b/mat/programčki/ničle.singlethreaded.c new file mode 100644 index 0000000..d2d1e0c --- /dev/null +++ b/mat/programčki/ničle.singlethreaded.c @@ -0,0 +1,74 @@ +#include <stdlib.h> // navdih za to je tedx predstavitev andreja bauerja z naslovom ničle +#include <stdio.h> // to je izris ničel littlewoodovih polinomov v manj kot 105 vrsticah +#include <gsl/gsl_poly.h> // prevod: gcc -Wall -Wextra -pedantic -lgsl -lm ničle.c +#include <sys/types.h> // program na poti "ime" izdela kvadratno PGM datoteko podane višine +#include <sys/stat.h> // za nadaljno obdelavo je uporabno, če je širina liha; ni pa nujno +#include <fcntl.h> // algoritem se da paralelizirati, ampak jaz imam samo en procesor +#include <sys/mman.h> // enotska krožnica je na četrtini podane širine +#include <string.h> +#include <unistd.h> +#include <gsl/gsl_errno.h> +void pripravi_koeficiente (double * izhod, int številka) { + while (številka) { + *izhod++ = številka & 1 ? 1 : -1; + številka >>= 1; + } +} +int main (int argc, char ** argv) { + int izven_slike = 0; + int nekonvergiranih = 0; + if (argc != 1+3) { + fprintf(stderr, "uporaba: %s stopnja ime širina\n", argv[0] ? argv[0] : "ničle"); + return 1; + } + int šir = atoi(argv[3]); + int fd; + if ((fd = open(argv[2], O_CREAT | O_RDWR, 00664)) == -1) { + perror("open"); + return 2; + } + if (ftruncate(fd, 128 + šir*šir) == -1) { + perror("ftruncate"); + if (close(fd) == -1) + perror("close"); + return 3; + } + void * p; + if ((p = mmap(NULL, 128 + šir*šir, PROT_READ|PROT_WRITE, MAP_SHARED, fd, 0)) == MAP_FAILED) { + perror("mmap"); + if (close(fd) == -1) + perror("close"); + return 4; + } + unsigned char * slika = (unsigned char *) p + 128; + memset(p, 0, 128 + šir*šir); + sprintf(p, "P5\n\n%59d\n%59d\n255\n", šir, šir); // precisely calculated with dc(1) (: + int stopnja = atoi(argv[1]); + double koeficienti[stopnja+1]; // kako prikladno! polinom nte stopnje ima n+1 členov + double ničle[2*stopnja]; // ima pa n ničer, 2n+0 so realni deli, 2n+1 pa imaginarni + gsl_set_error_handler_off(); + gsl_poly_complex_workspace * w = gsl_poly_complex_workspace_alloc(stopnja+1); + for (int i = 0; i < 1 << (stopnja+1); i++) { + pripravi_koeficiente(koeficienti, i); // noben člen ni 0, vsi so bodisi 1 bodisi -1 + if (gsl_poly_complex_solve(koeficienti, stopnja+1, w, ničle) != GSL_SUCCESS) + nekonvergiranih++; // uuu, lahko bi recimo narisali tiste, ki ne konver. + for (int j = 0; j < 2*stopnja; j += 2) { + int višina_na_sliki = šir/2 - ničle[j+1]*(šir/4); + int širina_na_sliki = šir/2 + ničle[j]*(šir/4); + if (višina_na_sliki > šir || širina_na_sliki > šir + || višina_na_sliki < 0 || širina_na_sliki < 0) { + izven_slike++; + continue; + } + slika[šir*višina_na_sliki+širina_na_sliki]++; + } + } + gsl_poly_complex_workspace_free(w); + if (munmap(p, 128 + šir*šir) == -1) + perror("munmap"); // nima smisla ukinjat programa, sistem si je sam kriv >:) + if (close(fd) == -1) + perror("close"); + printf("%d ničel je izven 2+2i (izven slike)\n%d polinomov ni konvergiralo\n", + izven_slike, nekonvergiranih); + return 0; +} |