Archiwum: Ogólna

Smogowy armageddon

Przez ostatnie dni media zdominował temat tragicznej jakości powietrza w polskich miastach. Nie jest to problem wirtualny, bo wystarczy wyjść na dwór by zobaczyć i poczuć marną jakość powietrza. Do wytworzenia smogu potrzebujemy dwóch czynników. Pierwszy to produkcja szkodliwych pyłów i zanieczyszczeń, która zimą jest mocno nasilona ze względu na sezon grzewczy, szczególnie w indywidualnych instalacjach C.O. z piecem na paliwo stało (drewno, węgiel), często opalanym odpadami. Elektrociepłownie, mimo, że opalane węglem wydalają proporcjonalnie mniej zanieczyszczeń na kW energii ze względu na zastosowane w instalacji filtry. Drugim, niezbędnym do powstania smogu warunkiem jest wystąpienie warstwy inwersyjnej w atmosferze, która blokuje wymianę powierza przy ziemi, z tym powyżej granicy inwersji. Jeżeli wysokość warstwy inwersyjnej wynosi 100-200 metrów, to znajdujemy się w zamkniętym “słoiku”, do którego wydalamy ogromne ilości zanieczyszczeń. Tyle skróconej teorii, przejdźmy do praktyki.

Wykorzystując bardzo duże stężenie pyłów w Warszawie postanowiłem sprawdzić jak problem wygląda z góry. Udało mi się wykonać dwa bardzo wymowne ujęcia:

Zdjęcie po lewej przedstawia kominy elektrociepłowni Siekierki i warszawskie drapacze chmur w tle. Zdjęcie po prawej to widok na centrum Warszawy z okolicy mostu Gdańskiego. Niestety zdjęcia nie pokazują chmur a właśnie zanieczyszczenia, którymi byliśmy zmuszeni oddychać.

Równolegle uruchomiłem pomiar stężenia pyłów w powietrzu na zewnątrz domu w Legionowie. Zakupiłem laserowy czujnik i połączyłem go z RaspberryPi. Czujnik wisi sobie bezpośrednio za oknem i komunikuje się z jeżyną przez interfejs szeregowy:

Każdy pomiar daje trzy wyniki odpowiadające pyłom o różnych frakcjach, w tym “najpopularniejszy” PM2.5. Program odbierający pomiary napisany w Pythonie:

import serial
from time import gmtime, strftime

port = serial.Serial("/dev/serial0", baudrate=9600, timeout=1.5)
data = port.read(32);

f = open("/home/pi/" + strftime("%Y-%m-%d", gmtime()) + ".txt", "a")
if ord(data[0]) == 66 and ord(data[1])==77:
	suma = 0
	for a in range(30):
		suma += ord(data[a])

	if suma == ord(data[30])*256+ord(data[31]):
		
		PM01 = str(ord(data[4])*256+ord(data[5]))
		PM25 = str(ord(data[6])*256+ord(data[7]))
		PM10 = str(ord(data[8])*256+ord(data[9]))
		str = strftime("%Y-%m-%d %H:%M:%S", gmtime()) + "\t" + PM01 + "\t" + PM25 + "\t" + PM10;
		f.write(str + "\n")
		print(str)
	
	
port.close()
f.close()

Powyższy kod uruchamiany jest co minutę. Dane zapisywane są w zwykłym pliku tekstowym:

2017-01-09 20:43:02	185	470	729
2017-01-09 20:44:02	181	440	682
2017-01-09 20:45:02	175	460	701
2017-01-09 20:46:02	172	445	725
2017-01-09 20:47:02	175	448	680
2017-01-09 20:48:02	177	449	725
2017-01-09 20:49:02	170	450	754
2017-01-09 20:50:02	176	458	746
2017-01-09 20:51:02	165	441	733
2017-01-09 20:52:02	174	434	707
2017-01-09 20:53:02	170	451	715
2017-01-09 20:54:02	170	436	702

Następnie raz na 15 minut uruchamiany jest kod rysujący wykres. Dla poprawy czytelności rysowane są średnie kroczące z pomiarów, a nie same pomiary:

import numpy
import ftplib
from ftplib import FTP
import sys
from matplotlib.dates import strpdate2num
from datetime import datetime, timedelta
import matplotlib

matplotlib.use('Agg')
import matplotlib.pyplot as plt
date = sys.argv[1];
X = numpy.loadtxt('/root/AIR/OLSZANKOWA50A/'+date+'.txt', delimiter='\t', dtype=object, converters={0: lambda x: datetime.strptime(x.decode("utf-8"), "%Y-%m-%d %H:%M:%S"), 1: numpy.float, 2: numpy.float, 3: numpy.float})

window = 900;
step = 60;
ts = min(X[:,0])
tss = datetime(ts.year, ts.month, ts.day)
t = tss;
te = t + timedelta(days=1)

T = []
PM01 = []
PM25 = []
PM10 = []
while t < te:
	t = t + timedelta(seconds=step)
	idx = numpy.logical_and(X[:,0] > (t - timedelta(seconds=window)), X[:,0] < (t + timedelta(seconds=window)))
	#print(t, sum(idx))
	T.append(t)
	if sum(idx) > 0:
		PM01.append(numpy.mean(X[idx,1]))
		PM25.append(numpy.mean(X[idx,2]))
		PM10.append(numpy.mean(X[idx,3]))
	else:
		PM01.append(numpy.nan)
		PM25.append(numpy.nan)
		PM10.append(numpy.nan)


x = []
label = []
for a in range(0,24,1):
	x.append(tss+timedelta(hours=a))
	if a%2 ==0:
		label.append((tss+timedelta(hours=a)).strftime('%H:%M:%S'))
	else:
		label.append('')

plt.figure(figsize=(20, 10))
plt.plot(T,PM01, 'b', lw=2, label='PM 1.00')
plt.plot(T,PM25, 'r', lw=4, label='PM 2.50')
plt.plot(T,PM10, 'g', lw=2, label='PM 10.0')
plt.xlabel('Time [UTC]')
plt.ylabel('ug/m3')
plt.legend(loc='upper left', shadow=True)
plt.xlim([ts, te])
plt.ylim([0, (numpy.ceil(max(PM10)/50)+1)*50])
plt.grid()
plt.title(date)
plt.xticks(x,label, rotation=0)
plt.savefig('/root/AIR/OLSZANKOWA50A/'+date+'.png',dpi=300)

File2Send = '/root/AIR/OLSZANKOWA50A/'+date+'.png'
Output_Directory = "/public_html/parking-legionowo.pl/AIR/" 


ftp = FTP('***************')
ftp.login('***************', '***************') 
file = open(File2Send, "rb") 
ftp.cwd(Output_Directory)
ftp.storbinary('STOR ' +date+'.png', file) 
print("STORing File now...")
ftp.quit() 
file.close() 
print("File transfered")

Po przygotowaniu wykresu jest on zapisywany do pliku png i wysyłany na serwer, gdzie można na bieżąco podglądać pomiary.

Tak wygląda normalny dzień “bez smogu” (te chwilowe wzrosty to prawdopodobnie momenty, gdy sąsiedzi dokładają węgla do starego pieca):

A tak wygląda dzień “smogowy” (akurat wtedy robiłem zdjęcia pokazane na początku):

Aktualne pomiary można zobaczyć na stronie sponsora: http://parking-legionowo.pl/AIR/

| Komentarze (2)

Nowy rok z lotu ptaka

Żeby zrobić ładnie zdjęcie dronem trzeba trafić na dobrą pogodę (nie może padać, zbyt mocno wiać, itp.) i dobrą porę dnia (najlepiej okolice wschodu i zachodu słońca, albo ładnie chmury). Są takie wydarzenia, które zdarzają się jeszcze rzadziej, a warte są uwiecznienia za pomocą drona: np. sylwester. Tylko raz w roku jest szansa sfotografowania i nagrania zmasowanego ataku fajerwerkowego. W tym roku w godzinę 0 panowała idealna do lotów pogoda: nie padało, nie wiało, a widoczność była rewelacyjna. Udało mi się wykonać kilka ładnych zdjęć i nagrać film:

Często pytacie mnie o fotografowanie nocą i długie ekspozycje za pomocą drona, więc postaram się niedługo odpowiedzieć Wam przez film na YouTube. Na razie zapraszam na film z nowego roku:

Jeżeli podobają się Wam moje zdjęcia i filmy to obserwujcie mnie na Instagramie (https://www.instagram.com/m.polkowski/) i subskrybujcie na YouTube (https://www.youtube.com/channel/UC0rcoI-FVLZjhIgTgfQpl9Q). Zapraszam!

| Komentarze

Świątecznie zadanie dla studentów

Co roku z okazji zbliżających się świąt proponuję studentom wstępu do programowania w pythonie na 1 roku w ramach kartkówki zadanie odmóżdżające polegające na przygotowaniu wykresu o tematyce świątecznej. Święta zbiegają się mniej więcej z wprowadzeniem do biblioteki matplotlib, więc i zadanie nie jest mocno naciągane od strony merytorycznej. Gdy dawałem to zadnie pierwszy raz byłem przekonany, że jest to z mojej strony prezent w postaci darmowego punktu do zaliczenia. Niestety co roku zdarzają się osoby, dla których wymyślenie współrzędnych kolejnych wierzchołków łamanej w kształcie choinki jest zbyt trudne.

Zadanie to ma jeszcze jeden ukryty cel – bardzo ładnie pokazuje czy dany student jest dobrym materiałem na fizyka-programistę:

  • Są osoby, które każdy element “rysunku” tworzą oddzielnie – nawet elementy choinki bywają narysowane jako pojedyncze proste. Gwiazdki / bombki również są w takich przypadkach rysowane przez wielokrotne wywołanie plot(). Autorzy takich prac z reguły mają duży problem z zaliczeniem tego prostego i wprowadzającego kursu.
  • Standardem są rozwiązania odrobinę bardziej przemyślane, gdzie wierzchołki łamanych i dodatkowe elementy są definiowane jako listy. O ile sam sposób rysowania może być przemyślany, o tyle współrzędnej dalej są wpisane ręcznie w kod programu. Większość autorów takich rozwiązań nie ma problemu z zaliczeniem przedmiotu na czwórkę.
  • Najrzadsze i najlepsze rozwiązania uwzględniają różnego rodzaju reguły matematyczne: np. funkcja generująca współrzędne łamanej tworzącej choinkę przyjmująca ilość gałęzi i wysokość jako parametry. Taka rozwiązanie spotkałem nie częściej niż raz na rok. Autorzy zaliczają zawsze na piątkę i nie mają problemów z pozostałymi (trudniejszymi) przedmiotami.

Najważniejsze w tym zadaniu jest nie sugerowania ani efektu końcowego, ani sposobu rozwiązania. Dzięki temu uruchamiana jest kreatywność (lub jej całkowity brak) i podświadome umiejętności programistyczne: czy dana osoba traktuje komputer jak ślepego wykonawce poleceń, czy jako maszynę która potrafi postępować na podstawie algorytmu.

Myślę, że mogłoby to być bardzo sprytne zadanie rekrutacyjne. Moi studenci dostawali na rozwiązanie 15-20 minut. Oto kilka przykładów:

| Komentarze

Na jakiej wysokości latają samoloty?

Analizując przepisy odnośnie ograniczeń lotów bezzałogowców w obrębie lotnisk zacząłem się zastanawiać na jakich faktycznie wysokościach latają samoloty rejsowe podczas startów i lądowań.

Wybrałem sobie jeden dzień (zeszły poniedziałek) i zebrałem informacje o wszystkich lotach nad Polską z tego dnia. Plik z danymi z całego dnia zajmuje jakieś 370 mb i wygląda mniej więcej tak:

bbc8f8f,3C5432,55.0600,11.8611,356,34000,354,5412,F-EKEB2,A306,D-AEAR,1480291202,LEJ,OSL,QY3314,0,0,BCS3314,0
bbc663a,394A0B,54.5979,11.8800,52,32000,401,0652,F-EKKL1,B77W,F-GSQL,1480291202,CDG,HND,AF274,0,0,AFR274,0
bbc879c,738076,49.0065,11.9041,120,37000,512,4753,F-ETHF1,B772,4X-ECF,1480291203,LHR,TLV,LY318,0,0,ELY318,0
bbc95e8,40066A,53.9157,12.4301,5,24000,349,6450,F-EDDT2,B752,G-BMRH,1480291202,LEJ,CPH,QY178,0,0,BCS178,0
bbbe339,,50.2147,12.6663,292,40025,385,3253,F-LKPD1,CL60,,1480291203,,,,0,0,CL60,0
bbc715e,406F01,49.8304,12.6811,16,35975,413,1000,F-LKVO1,A320,G-EZOZ,1480291203,BCN,SXF,U24536,0,64,EZY81TC,0
bbc86f4,424293,48.7195,13.1555,114,38975,516,2223,F-LKVO1,B788,VP-BBR,1480291203,LHR,GYD,J28,0,-64,AHY8,0
bbc858c,894070,50.5559,13.2681,112,35000,491,2246,F-LKVO1,A320,A9C-AO,1480291203,LHR,BAH,GF6,0,-64,GFA6,0

Mamy kod lotu, współrzędne geograficzne, wysokość, czas, lotnisko z którego wystartował, lotnisko docelowe, model samolotu itp. Komplet danych załadowałem do Matlaba. W pierwszej kolejności wyrysowałem sobie wszystkie dane:

mapa0

Sam ogrom informacji na mapie robi wrażenie, wyraźnie widać duży ruch w okolicach Okęcia. Mieliśmy rozmawiać o wysokościach, więc zaznaczam na mapie pozycje kolorami w zależności od wysokości lotu. Zielone to wysokości od 3 do 2 km, żółte od 2 do 1 km, pomarańczowe od 1 km do 500 m i czerwone poniżej 500 m:

mapa1

Widać, że obszar wokół lotnisk gdzie samoloty latają nisko jest stosunkowo duży. Dla każdego lotniska należałoby przeprowadzić oddzielną analizę. Ja zainteresowałem się Warszawskim Okęciem. Na wykresie naniosłem wysokość lotu (nad poziom terenu) w zależności od odległości w km od lotniska (dla Okęcia liczyłem odległość od przecięcia pasów startowych). Zbiór danych jest bardzo duży, a ich dokładność różna, ale wyraźnie widać różnicą pomiędzy profilem starów (czerwone) i lądowań (niebieskie):

wykres1

Wyraźnie widać, że w zależności od długości lotu samoloty osiągają różne wysokości przelotowe.

I to co najbardziej interesujące – zbliżenie na samo lotnisko:

wykres2

Wyraźnie widać, że podczas lądowań samoloty znajdują się dłużej i dalej na małej wysokości. Wydaje się to być zgodne ze zdrowym rozsądkiem. Jakby ktoś był zainteresowany mogę zrobić wykres dla innego lotniska. Bardzo ciekawe jest to zagięcie podczas lądowań równo 15 km od lotniska – może jakiś znawca tematu wyjaśni, czemu to tak wygląda. Wydawało by się, że powinno być łagodne rozpoczęcie zniżania a tu wydaje się gwałtowne.

Dla lubiących kod dziś mam fragment Matlaba rysujący wykres pokazany powyżej:

clear;clc;
T = readtable('all2.txt');

%%
airport1 = table2array(T(:,13));
airport2 = table2array(T(:,14));

flight_lon = table2array(T(:,4));
flight_lat = table2array(T(:,3));
flight_elev = table2array(T(:,6));
flight_code = table2array(T(:,1));
flight_elev = flight_elev * 0.3048;

%%
D = distance(52.165802, 20.967240, flight_lat(idx), flight_lon(idx));
D = deg2km(D);

idxs = strcmp(airport1,'WAW');
idxl = strcmp(airport2,'WAW');

Ds = deg2km(distance(52.165802, 20.967240, flight_lat(idxs), flight_lon(idxs)));
Dl = deg2km(distance(52.165802, 20.967240, flight_lat(idxl), flight_lon(idxl)));

clf
hold on
subplot(2,1,1)
plot(Dl,flight_elev(idxl)/1000,'b.')
title('Planes landing @ WAW')
xlabel('Distance [km]')
ylabel('Elevation [km]')
subplot(2,1,2)
plot(Ds,flight_elev(idxs)/1000,'r.')
title('Planes taking off @ WAW')
xlabel('Distance [km]')
ylabel('Elevation [km]')

print(gcf, '-r360', '-dpng', 'wykres1.png');

| Komentarze (2)

Prosty problem jako zadanie z programowania

Każdy pomysł na fajne zadanie z programowania jest bardzo cenny. Dobre zadanie dla studentów powinno być ciekawe i powiązane z jakimś życiowym problemem, tak aby łatwiej było się skupić na samej istocie programowania a nie na samym problemie. Sam chodziłem kiedyś na zajęcia, gdzie prowadzący wymyślał abstrakcyjne problemy, które łatwo było oprogramować, ale bardzo trudno było zrozumieć o co chodzi i po co to wszystko. Wczoraj otworzyłem nowy wielopak zapałek i okazały się to być te z klasyczną łamigłówką z przekładaniem jednej zapałki:

img_1110

Pomyślałem, że mogło by to być fajne zadanie dla studentów I roku. Ja zazwyczaj uczę Pythona, więc podjąłem próbę napisania programu. W pierwszej wersji zadaniem programu jest po prostu znalezienie rozwiązania (lub rozwiązań) pojedynczego problemu:


plus = dict()
plus[0] = [8]
plus[1] = [7]
plus[2] = []
plus[3] = [9]
plus[4] = []
plus[5] = [9]
plus[6] = [8]
plus[7] = []
plus[8] = []
plus[9] = [8]

minus = dict()
minus[0] = []
minus[1] = []
minus[2] = []
minus[3] = []
minus[4] = []
minus[5] = []
minus[6] = []
minus[7] = [1]
minus[8] = [0, 6, 9]
minus[9] = [3, 5]

swap = dict()
swap[0] = [6, 9]
swap[1] = []
swap[2] = [3]
swap[3] = [2, 5]
swap[4] = []
swap[5] = [3]
swap[6] = [0, 9]
swap[7] = []
swap[8] = []
swap[9] = [0, 6]


def test_single(a, o1, b, o2, c):
	if o2 == "=":
		if o1 == "+":
			return a+b==c
		elif o1 == "-":
			return a-b==c
		else:
			print('ERROR')
	elif o1 == "=":
		if o2 == "+":
			return a+b==c
		elif o2 == "-":
			return a-b==c
		else:
			print('ERROR')
	else:
		print('ERROR')
def print_result(a, o1, b, o2, c,debug):
	print("{5}: {0}{1}{2}{3}{4}:".format(a,o1,b,o2,c,debug), test_single(a, o1, b, o2, c));


def all_comb(a, o1, b, o2, c):
	#print("SWAPY")
	for a1 in swap[a]:
		print_result(a1, o1, b, o2, c,1)
	for b1 in swap[b]:
		print_result(a, o1, b1, o2, c,2)
	for c1 in swap[c]:
		print_result(a, o1, b, o2, c1,3)
	if o1 == '-':
		print_result(a, "=", b, "-", c,4)
	#print("A minus")
	for a1 in minus[a]:
		for b1 in plus[b]:
			print_result(a1, o1, b1, o2, c,5)
		for c1 in plus[c]:
			print_result(a1, o1, b, o2, c1,6)
		if o1 == "-":
			print_result(a1, "+", b, o2, c,7)
		if o2 == "-":
			print_result(a1, o1, b, "+", c,8)
	#print("B minus")
	for b1 in minus[b]:
		for a1 in plus[a]:
			print_result(a1, o1, b1, o2, c,9)
		for c1 in plus[c]:
			print_result(a, o1, b1, o2, c1,10)
		if o1 == "-":
			print_result(a, "+", b1, o2, c,11)
		if o2 == "-":
			print_result(a, o1, b1, "+", c,12)
	#print("C minus")
	for c1 in minus[c]:
		for b1 in plus[b]:
			print_result(a, o1, b1, o2, c1,13)
		for a1 in plus[a]:
			print_result(a1, o1, b, o2, c1,14)
		if o1 == "-":
			print_result(a, "+", b, o2, c1,15)
		if o2 == "-":
			print_result(a, o1, b, "+", c1,16)

	#print("+ minus")
	if o1 == "+":
		for a1 in plus[a]:
			print_result(a1, "-", b, o2, c,17)
		for b1 in plus[b]:
			print_result(a, "-", b1, o2, c,18)
		for c1 in plus[c]:
			print_result(a, "-", b, o2, c1,19)


all_comb(1,"+",2,"=",8)

Nie jest to rozwiązanie w żaden sposób zoptymalizowane. Byłbym zadowolony, gdyby na zajęciach udało się dojść do takiego rozwiązania. Dla równania 1+2=8 program zwraca taki wynik:

screenshot_43

Działanie programu opiera się po pierwsze na 3 zdefiniowanych słownikach, które określają na jaką cyfrę można zmienić daną cyfrę po odjęciu, dodaniu i przełożeniu jednej zapałki. Po drugie program przewiduje wszystkie przypadki zamian: np. odjęcie zapałki z pierwszej cyfry i dodanie do cyfry drugiej lub trzeciej. W sumie przewidziane jest 19 przypadków (może da się to zoptymalizować?).

Drugą, równie ciekawą wersją programu jest wygenerowanie wszystkich możliwych łamigłówek. Tu należy zdecydować, czy rozważamy te równania, które są poprawne na wejściu. Ja zdecydowałem ich nie pomijać, bo wymagamy poprawności po przełożeniu jednej zapałki. Kod wygląda tak:


plus = dict()
plus[0] = [8]
plus[1] = [7]
plus[2] = []
plus[3] = [9]
plus[4] = []
plus[5] = [9]
plus[6] = [8]
plus[7] = []
plus[8] = []
plus[9] = [8]

minus = dict()
minus[0] = []
minus[1] = []
minus[2] = []
minus[3] = []
minus[4] = []
minus[5] = []
minus[6] = []
minus[7] = [1]
minus[8] = [0, 6, 9]
minus[9] = [3, 5]

swap = dict()
swap[0] = [6, 9]
swap[1] = []
swap[2] = [3]
swap[3] = [2, 5]
swap[4] = []
swap[5] = [3]
swap[6] = [0, 9]
swap[7] = []
swap[8] = []
swap[9] = [0, 6]


def test_single(a, o1, b, o2, c):
	if o2 == "=":
		if o1 == "+":
			return a+b==c
		elif o1 == "-":
			return a-b==c
		else:
			print('ERROR')
	elif o1 == "=":
		if o2 == "+":
			return a+b==c
		elif o2 == "-":
			return a-b==c
		else:
			print('ERROR')
	else:
		print('ERROR')
def print_result(wa,wo1,wb,wo2,wc,a, o1, b, o2, c):
	if test_single(a, o1, b, o2, c):
		print("{0}{1}{2}{3}{4} -> {5}{6}{7}{8}{9}".format(wa,wo1,wb,wo2,wc,a,o1,b,o2,c))
		return 1
	else:
		return 0;


def all_comb(a, o1, b, o2, c):
	count = 0
	for a1 in swap[a]:
		count = count + print_result(a,o1,b,o2,c,a1, o1, b, o2, c)
	for b1 in swap[b]:
		count = count + print_result(a,o1,b,o2,c,a, o1, b1, o2, c)
	for c1 in swap[c]:
		count = count + print_result(a,o1,b,o2,c,a, o1, b, o2, c1)
	if o1 == '-':
		count = count + print_result(a,o1,b,o2,c,a, "=", b, "-", c)
	for a1 in minus[a]:
		for b1 in plus[b]:
			count = count + print_result(a,o1,b,o2,c,a1, o1, b1, o2, c)
		for c1 in plus[c]:
			count = count + print_result(a,o1,b,o2,c,a1, o1, b, o2, c1)
		if o1 == "-":
			count = count + print_result(a,o1,b,o2,c,a1, "+", b, o2, c)
		if o2 == "-":
			count = count + print_result(a,o1,b,o2,c,a1, o1, b, "+", c)
	for b1 in minus[b]:
		for a1 in plus[a]:
			count = count + print_result(a,o1,b,o2,c,a1, o1, b1, o2, c)
		for c1 in plus[c]:
			count = count + print_result(a,o1,b,o2,c,a, o1, b1, o2, c1)
		if o1 == "-":
			count = count + print_result(a,o1,b,o2,c,a, "+", b1, o2, c)
		if o2 == "-":
			count = count + print_result(a,o1,b,o2,c,a, o1, b1, "+", c)
	for c1 in minus[c]:
		for b1 in plus[b]:
			count = count + print_result(a,o1,b,o2,c,a, o1, b1, o2, c1)
		for a1 in plus[a]:
			count = count + print_result(a,o1,b,o2,c,a1, o1, b, o2, c1)
		if o1 == "-":
			count = count + print_result(a,o1,b,o2,c,a, "+", b, o2, c1)
		if o2 == "-":
			count = count + print_result(a,o1,b,o2,c,a, o1, b, "+", c1)
	if o1 == "+":
		for a1 in plus[a]:
			count = count + print_result(a,o1,b,o2,c,a1, "-", b, o2, c)
		for b1 in plus[b]:
			count = count + print_result(a,o1,b,o2,c,a, "-", b1, o2, c)
		for c1 in plus[c]:
			count = count + print_result(a,o1,b,o2,c,a, "-", b, o2, c1)
	return count

for to1 in ["+", "-"]:
	for ta in range(10):
		for tb in range(10):
			for tc in range(10):
				tcount = all_comb(ta,to1,tb,"=",tc)
				#if tcount > 2:
				#	print("{0}{1}{2}{3}{4}".format(ta,to1,tb,"=",tc), tcount)

Wynik jest obliczany w mgnieniu oka. Dla rozwiązania dla łamigłówek z obrazka są następujące:

screenshot_44

Lista wszystkich kombinacji dla potomnych:

0+0=6 -> 6+0=6
0+0=6 -> 0+6=6
0+0=6 -> 0+0=0
0+0=8 -> 8-0=8
0+0=9 -> 9+0=9
0+0=9 -> 0+9=9
0+0=9 -> 0+0=0
0+1=7 -> 6+1=7
0+1=7 -> 8-1=7
0+1=8 -> 8+1=9
0+2=3 -> 0+3=3
0+2=3 -> 0+2=2
0+2=6 -> 8-2=6
0+2=8 -> 6+2=8
0+3=2 -> 0+2=2
0+3=2 -> 0+3=3
0+3=5 -> 0+5=5
0+3=5 -> 0+3=3
0+3=5 -> 8-3=5
0+3=8 -> 0+9=9
0+3=9 -> 6+3=9
0+4=4 -> 8-4=4
0+5=3 -> 0+3=3
0+5=3 -> 0+5=5
0+5=3 -> 8-5=3
0+5=8 -> 0+9=9
0+6=0 -> 0+0=0
0+6=0 -> 0+6=6
0+6=2 -> 8-6=2
0+6=9 -> 0+9=9
0+6=9 -> 0+6=6
0+7=1 -> 8-7=1
0+7=9 -> 8+1=9
0+8=0 -> 8-8=0
0+8=3 -> 0+9=9
0+8=5 -> 0+9=9
0+8=8 -> 8+0=8
0+9=0 -> 0+0=0
0+9=0 -> 0+9=9
0+9=6 -> 0+6=6
0+9=6 -> 0+9=9
1+0=7 -> 1+6=7
1+0=7 -> 7-0=7
1+0=8 -> 1+8=9
1+1=3 -> 1+1=2
1+1=6 -> 7-1=6
1+2=2 -> 1+2=3
1+2=4 -> 1+3=4
1+2=5 -> 1+2=3
1+2=5 -> 7-2=5
1+2=8 -> 7+2=9
1+3=3 -> 1+2=3
1+3=4 -> 7-3=4
1+3=6 -> 1+5=6
1+4=3 -> 1+4=5
1+4=3 -> 7-4=3
1+5=0 -> 1+5=6
1+5=2 -> 7-5=2
1+5=4 -> 1+3=4
1+5=9 -> 1+5=6
1+6=1 -> 1+0=1
1+6=1 -> 7-6=1
1+6=8 -> 1+8=9
1+7=0 -> 7-7=0
1+7=8 -> 7+1=8
1+8=0 -> 1+8=9
1+8=1 -> 1+6=7
1+8=6 -> 1+8=9
1+8=7 -> 7+0=7
1+9=1 -> 1+0=1
1+9=7 -> 1+6=7
1+9=8 -> 1+8=9
2+0=3 -> 3+0=3
2+0=3 -> 2+0=2
2+0=8 -> 2+6=8
2+1=2 -> 2+1=3
2+1=4 -> 3+1=4
2+1=5 -> 2+1=3
2+1=8 -> 2+7=9
2+2=5 -> 3+2=5
2+2=5 -> 2+3=5
2+3=3 -> 2+3=5
2+3=4 -> 2+2=4
2+3=6 -> 3+3=6
2+3=7 -> 2+5=7
2+4=0 -> 2+4=6
2+4=7 -> 3+4=7
2+4=9 -> 2+4=6
2+5=5 -> 2+3=5
2+5=8 -> 3+5=8
2+6=2 -> 2+0=2
2+6=9 -> 3+6=9
2+7=0 -> 2+7=9
2+7=6 -> 2+7=9
2+8=0 -> 2+6=8
2+8=6 -> 2+6=8
2+8=9 -> 2+6=8
2+9=1 -> 2+5=7
2+9=2 -> 2+0=2
2+9=8 -> 2+6=8
3+0=2 -> 2+0=2
3+0=2 -> 3+0=3
3+0=5 -> 5+0=5
3+0=5 -> 3+0=3
3+0=8 -> 9+0=9
3+0=9 -> 3+6=9
3+0=9 -> 9-0=9
3+1=3 -> 2+1=3
3+1=6 -> 5+1=6
3+1=8 -> 9-1=8
3+2=3 -> 3+2=5
3+2=4 -> 2+2=4
3+2=6 -> 3+3=6
3+2=7 -> 5+2=7
3+2=7 -> 9-2=7
3+3=0 -> 3+3=6
3+3=5 -> 2+3=5
3+3=5 -> 3+2=5
3+3=6 -> 9-3=6
3+3=8 -> 5+3=8
3+3=8 -> 3+5=8
3+3=9 -> 3+3=6
3+4=5 -> 9-4=5
3+4=6 -> 2+4=6
3+4=9 -> 5+4=9
3+5=4 -> 9-5=4
3+5=6 -> 3+3=6
3+5=7 -> 2+5=7
3+6=0 -> 3+6=9
3+6=3 -> 3+0=3
3+6=3 -> 9-6=3
3+6=6 -> 3+6=9
3+6=8 -> 2+6=8
3+7=2 -> 9-7=2
3+7=9 -> 2+7=9
3+8=1 -> 9-8=1
3+8=3 -> 3+6=9
3+8=5 -> 3+6=9
3+8=9 -> 9+0=9
3+9=0 -> 3+5=8
3+9=0 -> 9-9=0
3+9=3 -> 3+0=3
3+9=6 -> 3+5=8
3+9=9 -> 3+6=9
3+9=9 -> 3+5=8
4+1=3 -> 4+1=5
4+2=0 -> 4+2=6
4+2=7 -> 4+3=7
4+2=9 -> 4+2=6
4+3=6 -> 4+2=6
4+3=9 -> 4+5=9
4+5=0 -> 4+5=9
4+5=6 -> 4+5=9
4+5=7 -> 4+3=7
4+6=4 -> 4+0=4
4+9=1 -> 4+3=7
4+9=3 -> 4+5=9
4+9=4 -> 4+0=4
4+9=5 -> 4+5=9
5+0=3 -> 3+0=3
5+0=3 -> 5+0=5
5+0=8 -> 9+0=9
5+0=9 -> 9-0=9
5+1=0 -> 5+1=6
5+1=4 -> 3+1=4
5+1=8 -> 9-1=8
5+1=9 -> 5+1=6
5+2=5 -> 3+2=5
5+2=7 -> 9-2=7
5+2=8 -> 5+3=8
5+3=6 -> 3+3=6
5+3=6 -> 9-3=6
5+3=7 -> 5+2=7
5+4=0 -> 5+4=9
5+4=5 -> 9-4=5
5+4=6 -> 5+4=9
5+4=7 -> 3+4=7
5+5=4 -> 9-5=4
5+5=8 -> 3+5=8
5+5=8 -> 5+3=8
5+6=3 -> 9-6=3
5+6=5 -> 5+0=5
5+6=9 -> 3+6=9
5+7=2 -> 9-7=2
5+8=1 -> 9-8=1
5+8=9 -> 9+0=9
5+9=0 -> 5+3=8
5+9=0 -> 9-9=0
5+9=5 -> 5+0=5
5+9=6 -> 5+3=8
5+9=9 -> 5+3=8
6+0=0 -> 0+0=0
6+0=0 -> 6+0=6
6+0=8 -> 8-0=8
6+0=9 -> 9+0=9
6+0=9 -> 6+0=6
6+1=1 -> 0+1=1
6+1=7 -> 8-1=7
6+1=8 -> 8+1=9
6+2=2 -> 0+2=2
6+2=6 -> 8-2=6
6+2=9 -> 6+3=9
6+3=0 -> 6+3=9
6+3=3 -> 0+3=3
6+3=5 -> 8-3=5
6+3=6 -> 6+3=9
6+3=8 -> 6+2=8
6+4=4 -> 0+4=4
6+4=4 -> 8-4=4
6+5=3 -> 8-5=3
6+5=5 -> 0+5=5
6+5=9 -> 6+3=9
6+6=2 -> 8-6=2
6+6=6 -> 0+6=6
6+6=6 -> 6+0=6
6+7=1 -> 6+1=7
6+7=1 -> 8-7=1
6+7=7 -> 0+7=7
6+7=9 -> 8+1=9
6+8=0 -> 8-8=0
6+8=8 -> 0+8=8
6+8=8 -> 8+0=8
6+9=3 -> 6+3=9
6+9=5 -> 6+3=9
6+9=6 -> 6+0=6
6+9=9 -> 0+9=9
7+0=1 -> 7-0=7
7+0=9 -> 1+8=9
7+1=0 -> 7-7=0
7+1=8 -> 1+7=8
7+2=0 -> 7+2=9
7+2=6 -> 7+2=9
7+3=9 -> 7+2=9
7+6=1 -> 1+6=7
7+6=7 -> 7+0=7
7+6=9 -> 1+8=9
7+7=0 -> 1+7=8
7+7=0 -> 7+1=8
7+7=6 -> 1+7=8
7+7=6 -> 7+1=8
7+7=9 -> 1+7=8
7+7=9 -> 7+1=8
7+8=1 -> 7+0=7
7+8=3 -> 1+8=9
7+8=5 -> 1+8=9
7+9=7 -> 7+0=7
7+9=9 -> 1+8=9
8+0=0 -> 8-8=0
8+0=0 -> 8-0=8
8+0=3 -> 9+0=9
8+0=5 -> 9+0=9
8+0=6 -> 8-0=8
8+0=8 -> 0+8=8
8+0=9 -> 8-0=8
8+1=0 -> 8+1=9
8+1=1 -> 6+1=7
8+1=1 -> 8-7=1
8+1=1 -> 8-1=7
8+1=6 -> 8+1=9
8+1=7 -> 0+7=7
8+2=0 -> 6+2=8
8+2=6 -> 6+2=8
8+2=9 -> 6+2=8
8+3=3 -> 6+3=9
8+3=5 -> 6+3=9
8+3=9 -> 0+9=9
8+5=9 -> 0+9=9
8+6=0 -> 8-8=0
8+6=8 -> 8+0=8
8+6=8 -> 0+8=8
8+7=1 -> 0+7=7
8+7=3 -> 8+1=9
8+7=5 -> 8+1=9
8+8=0 -> 0+8=8
8+8=0 -> 8+0=8
8+8=6 -> 0+8=8
8+8=6 -> 8+0=8
8+8=9 -> 0+8=8
8+8=9 -> 8+0=8
8+9=0 -> 8-8=0
8+9=3 -> 0+9=9
8+9=5 -> 0+9=9
8+9=8 -> 8+0=8
8+9=8 -> 0+8=8
9+0=0 -> 0+0=0
9+0=0 -> 9+0=9
9+0=1 -> 9-8=1
9+0=3 -> 9-0=9
9+0=5 -> 9-0=9
9+0=6 -> 6+0=6
9+0=6 -> 9+0=9
9+0=8 -> 8-0=8
9+1=0 -> 9-1=8
9+1=1 -> 0+1=1
9+1=2 -> 9-7=2
9+1=6 -> 9-1=8
9+1=7 -> 6+1=7
9+1=7 -> 8-1=7
9+1=8 -> 8+1=9
9+1=9 -> 9-1=8
9+2=1 -> 5+2=7
9+2=1 -> 9-2=7
9+2=2 -> 0+2=2
9+2=6 -> 8-2=6
9+2=8 -> 6+2=8
9+3=0 -> 5+3=8
9+3=0 -> 9-9=0
9+3=3 -> 0+3=3
9+3=5 -> 8-3=5
9+3=6 -> 5+3=8
9+3=9 -> 6+3=9
9+3=9 -> 5+3=8
9+4=1 -> 3+4=7
9+4=3 -> 5+4=9
9+4=4 -> 0+4=4
9+4=4 -> 8-4=4
9+4=5 -> 5+4=9
9+5=0 -> 3+5=8
9+5=0 -> 9-9=0
9+5=3 -> 8-5=3
9+5=5 -> 0+5=5
9+5=6 -> 3+5=8
9+5=9 -> 3+5=8
9+6=1 -> 9-8=1
9+6=2 -> 8-6=2
9+6=3 -> 3+6=9
9+6=5 -> 3+6=9
9+6=6 -> 0+6=6
9+6=9 -> 9+0=9
9+7=1 -> 8-7=1
9+7=7 -> 0+7=7
9+7=9 -> 8+1=9
9+8=0 -> 8-8=0
9+8=3 -> 9+0=9
9+8=5 -> 9+0=9
9+8=8 -> 0+8=8
9+8=8 -> 8+0=8
9+9=1 -> 9-8=1
9+9=9 -> 0+9=9
9+9=9 -> 9+0=9
0-0=0 -> 0=0-0
0-0=6 -> 6-0=6
0-0=6 -> 0-0=0
0-0=8 -> 0+0=0
0-0=9 -> 9-0=9
0-0=9 -> 0-0=0
0-1=5 -> 6-1=5
0-1=7 -> 0+1=1
0-1=8 -> 9-1=8
0-2=4 -> 6-2=4
0-2=7 -> 9-2=7
0-2=8 -> 8-2=6
0-3=3 -> 6-3=3
0-3=6 -> 9-3=6
0-3=9 -> 0+3=3
0-3=9 -> 8-3=5
0-4=2 -> 6-4=2
0-4=5 -> 9-4=5
0-5=1 -> 6-5=1
0-5=4 -> 9-5=4
0-5=9 -> 8-5=3
0-5=9 -> 0+5=5
0-6=0 -> 6-6=0
0-6=0 -> 0-0=0
0-6=3 -> 9-6=3
0-6=8 -> 0+6=6
0-7=1 -> 0+1=1
0-7=2 -> 9-7=2
0-7=7 -> 8-1=7
0-7=7 -> 8-7=1
0-8=0 -> 0+0=0
0-8=1 -> 9-8=1
0-8=2 -> 8-6=2
0-8=6 -> 0+6=6
0-8=8 -> 8-0=8
0-8=8 -> 8-8=0
0-8=9 -> 0+9=9
0-9=0 -> 9-9=0
0-9=0 -> 0-0=0
0-9=3 -> 0+3=3
0-9=3 -> 8-5=3
0-9=5 -> 8-3=5
0-9=5 -> 0+5=5
0-9=8 -> 0+9=9
1-0=1 -> 1=0-1
1-0=7 -> 1+0=1
1-1=0 -> 1=1-0
1-1=6 -> 1-1=0
1-1=8 -> 7-1=6
1-1=9 -> 1-1=0
1-2=9 -> 1+2=3
1-2=9 -> 7-2=5
1-4=9 -> 7-4=3
1-4=9 -> 1+4=5
1-5=8 -> 1+5=6
1-6=1 -> 1-0=1
1-6=7 -> 7-6=1
1-7=2 -> 1+1=2
1-7=6 -> 7-1=6
1-7=8 -> 7-7=0
1-8=1 -> 1+0=1
1-8=1 -> 7-6=1
1-8=7 -> 7-0=7
1-8=7 -> 1+6=7
1-8=8 -> 1+8=9
1-9=1 -> 1-0=1
1-9=2 -> 7-5=2
1-9=4 -> 7-3=4
1-9=4 -> 1+3=4
1-9=6 -> 1+5=6
2-0=2 -> 2=0-2
2-0=3 -> 3-0=3
2-0=3 -> 2-0=2
2-1=1 -> 2=1-1
2-1=2 -> 3-1=2
2-1=9 -> 2+1=3
2-2=0 -> 2=2-0
2-2=1 -> 3-2=1
2-2=6 -> 2-2=0
2-2=9 -> 2-2=0
2-3=0 -> 3-3=0
2-3=0 -> 2-2=0
2-3=9 -> 2+3=5
2-4=8 -> 2+4=6
2-6=2 -> 2-0=2
2-7=3 -> 2+1=3
2-7=8 -> 2+7=9
2-8=2 -> 2+0=2
2-8=8 -> 2+6=8
2-9=2 -> 2-0=2
2-9=5 -> 2+3=5
2-9=7 -> 2+5=7
3-0=2 -> 2-0=2
3-0=2 -> 3-0=3
3-0=3 -> 3=0-3
3-0=5 -> 5-0=5
3-0=5 -> 3-0=3
3-0=8 -> 9-0=9
3-0=9 -> 3+0=3
3-1=1 -> 2-1=1
3-1=2 -> 3=1-2
3-1=3 -> 3-1=2
3-1=4 -> 5-1=4
3-2=0 -> 2-2=0
3-2=0 -> 3-3=0
3-2=1 -> 3=2-1
3-2=3 -> 5-2=3
3-2=9 -> 3+2=5
3-3=0 -> 3=3-0
3-3=1 -> 3-2=1
3-3=2 -> 5-3=2
3-3=6 -> 3-3=0
3-3=8 -> 9-3=6
3-3=8 -> 3+3=6
3-3=9 -> 3-3=0
3-4=1 -> 5-4=1
3-4=9 -> 9-4=5
3-5=0 -> 5-5=0
3-5=0 -> 3-3=0
3-6=3 -> 3-0=3
3-6=8 -> 3+6=9
3-6=9 -> 9-6=3
3-7=4 -> 3+1=4
3-7=8 -> 9-1=8
3-8=0 -> 9-9=0
3-8=3 -> 3+0=3
3-8=3 -> 9-6=3
3-8=7 -> 9-8=1
3-8=9 -> 9-0=9
3-8=9 -> 3+6=9
3-9=3 -> 3-0=3
3-9=4 -> 9-5=4
3-9=6 -> 9-3=6
3-9=6 -> 3+3=6
3-9=8 -> 3+5=8
3-9=8 -> 9-9=0
4-0=4 -> 4=0-4
4-1=2 -> 4-1=3
4-1=3 -> 4=1-3
4-1=5 -> 4-1=3
4-1=9 -> 4+1=5
4-2=1 -> 4-3=1
4-2=2 -> 4=2-2
4-2=3 -> 4-2=2
4-2=8 -> 4+2=6
4-3=1 -> 4=3-1
4-3=2 -> 4-2=2
4-4=0 -> 4=4-0
4-4=6 -> 4-4=0
4-4=9 -> 4-4=0
4-5=1 -> 4-3=1
4-5=8 -> 4+5=9
4-6=4 -> 4-0=4
4-7=5 -> 4+1=5
4-8=4 -> 4+0=4
4-9=4 -> 4-0=4
4-9=7 -> 4+3=7
4-9=9 -> 4+5=9
5-0=3 -> 3-0=3
5-0=3 -> 5-0=5
5-0=5 -> 5=0-5
5-0=8 -> 9-0=9
5-0=9 -> 5+0=5
5-1=2 -> 3-1=2
5-1=4 -> 5=1-4
5-1=8 -> 5+1=6
5-2=1 -> 3-2=1
5-2=2 -> 5-3=2
5-2=2 -> 5-2=3
5-2=3 -> 5=2-3
5-2=5 -> 5-2=3
5-3=0 -> 3-3=0
5-3=0 -> 5-5=0
5-3=2 -> 5=3-2
5-3=3 -> 5-2=3
5-3=3 -> 5-3=2
5-3=8 -> 9-3=6
5-4=1 -> 5=4-1
5-4=8 -> 5+4=9
5-4=9 -> 9-4=5
5-5=0 -> 5=5-0
5-5=2 -> 5-3=2
5-5=6 -> 5-5=0
5-5=9 -> 5-5=0
5-6=5 -> 5-0=5
5-6=9 -> 9-6=3
5-7=6 -> 5+1=6
5-7=8 -> 9-1=8
5-8=0 -> 9-9=0
5-8=3 -> 9-6=3
5-8=5 -> 5+0=5
5-8=7 -> 9-8=1
5-8=9 -> 9-0=9
5-9=4 -> 9-5=4
5-9=5 -> 5-0=5
5-9=6 -> 9-3=6
5-9=8 -> 5+3=8
5-9=8 -> 9-9=0
6-0=0 -> 0-0=0
6-0=0 -> 6-6=0
6-0=0 -> 6-0=6
6-0=6 -> 6=0-6
6-0=8 -> 6+0=6
6-0=9 -> 9-0=9
6-0=9 -> 6-0=6
6-1=3 -> 6-1=5
6-1=5 -> 6=1-5
6-1=8 -> 9-1=8
6-2=3 -> 6-3=3
6-2=4 -> 6=2-4
6-2=7 -> 9-2=7
6-2=8 -> 8-2=6
6-3=1 -> 6-5=1
6-3=2 -> 6-3=3
6-3=3 -> 6=3-3
6-3=4 -> 6-2=4
6-3=5 -> 6-3=3
6-3=6 -> 9-3=6
6-3=8 -> 6+3=9
6-3=9 -> 8-3=5
6-4=2 -> 6=4-2
6-4=3 -> 6-4=2
6-4=5 -> 9-4=5
6-5=1 -> 6=5-1
6-5=3 -> 6-3=3
6-5=4 -> 9-5=4
6-5=9 -> 8-5=3
6-6=0 -> 6=6-0
6-6=3 -> 9-6=3
6-6=6 -> 6-0=6
6-6=6 -> 6-6=0
6-6=9 -> 6-6=0
6-7=2 -> 9-7=2
6-7=7 -> 8-1=7
6-7=7 -> 6+1=7
6-7=7 -> 8-7=1
6-8=1 -> 9-8=1
6-8=2 -> 8-6=2
6-8=6 -> 6+0=6
6-8=8 -> 8-0=8
6-8=8 -> 8-8=0
6-9=0 -> 9-9=0
6-9=0 -> 6-6=0
6-9=3 -> 8-5=3
6-9=5 -> 8-3=5
6-9=6 -> 6-0=6
6-9=9 -> 6+3=9
7-0=1 -> 7-6=1
7-0=1 -> 1+0=1
7-0=7 -> 7=0-7
7-1=0 -> 7-1=6
7-1=2 -> 1+1=2
7-1=6 -> 7=1-6
7-1=8 -> 7-7=0
7-1=9 -> 7-1=6
7-2=3 -> 7-2=5
7-2=3 -> 1+2=3
7-2=4 -> 7-3=4
7-2=5 -> 7=2-5
7-2=8 -> 7+2=9
7-3=2 -> 7-5=2
7-3=4 -> 7=3-4
7-3=4 -> 1+3=4
7-3=5 -> 7-2=5
7-4=2 -> 7-4=3
7-4=3 -> 7=4-3
7-4=5 -> 7-4=3
7-4=5 -> 1+4=5
7-5=2 -> 7=5-2
7-5=3 -> 7-5=2
7-5=4 -> 7-3=4
7-5=6 -> 1+5=6
7-6=1 -> 7=6-1
7-6=7 -> 7-0=7
7-6=7 -> 1+6=7
7-7=0 -> 7=7-0
7-7=6 -> 7-7=0
7-7=8 -> 1+7=8
7-7=8 -> 7+1=8
7-7=9 -> 7-7=0
7-8=1 -> 7-0=7
7-8=7 -> 7+0=7
7-8=9 -> 1+8=9
7-9=1 -> 7-6=1
7-9=7 -> 7-0=7
8-0=0 -> 0+0=0
8-0=1 -> 9-8=1
8-0=2 -> 8-6=2
8-0=3 -> 9-0=9
8-0=5 -> 9-0=9
8-0=6 -> 6+0=6
8-0=8 -> 8=0-8
8-0=8 -> 8-8=0
8-0=9 -> 9+0=9
8-1=0 -> 9-1=8
8-1=1 -> 0+1=1
8-1=2 -> 9-7=2
8-1=6 -> 9-1=8
8-1=7 -> 8=1-7
8-1=7 -> 6+1=7
8-1=7 -> 8-7=1
8-1=8 -> 8+1=9
8-1=9 -> 9-1=8
8-2=0 -> 8-2=6
8-2=1 -> 9-2=7
8-2=2 -> 0+2=2
8-2=5 -> 8-3=5
8-2=6 -> 8=2-6
8-2=8 -> 6+2=8
8-2=9 -> 8-2=6
8-3=0 -> 9-9=0
8-3=3 -> 8-5=3
8-3=3 -> 8-3=5
8-3=3 -> 0+3=3
8-3=5 -> 8=3-5
8-3=6 -> 8-2=6
8-3=9 -> 6+3=9
8-4=4 -> 8=4-4
8-4=4 -> 0+4=4
8-5=0 -> 9-9=0
8-5=2 -> 8-5=3
8-5=3 -> 8=5-3
8-5=5 -> 8-3=5
8-5=5 -> 8-5=3
8-5=5 -> 0+5=5
8-6=1 -> 9-8=1
8-6=2 -> 8=6-2
8-6=3 -> 8-6=2
8-6=6 -> 0+6=6
8-6=8 -> 8-0=8
8-6=8 -> 8-8=0
8-7=1 -> 8=7-1
8-7=1 -> 8-1=7
8-7=7 -> 0+7=7
8-7=9 -> 8+1=9
8-8=0 -> 8=8-0
8-8=0 -> 8-0=8
8-8=6 -> 8-8=0
8-8=6 -> 8-0=8
8-8=8 -> 0+8=8
8-8=8 -> 8+0=8
8-8=9 -> 8-8=0
8-8=9 -> 8-0=8
8-9=1 -> 9-8=1
8-9=2 -> 8-6=2
8-9=8 -> 8-0=8
8-9=8 -> 8-8=0
8-9=9 -> 0+9=9
9-0=0 -> 0-0=0
9-0=0 -> 9-9=0
9-0=0 -> 9-0=9
9-0=3 -> 9-6=3
9-0=3 -> 3+0=3
9-0=5 -> 5+0=5
9-0=6 -> 6-0=6
9-0=6 -> 9-0=9
9-0=7 -> 9-8=1
9-0=8 -> 9+0=9
9-0=9 -> 9=0-9
9-1=4 -> 3+1=4
9-1=5 -> 6-1=5
9-1=6 -> 5+1=6
9-1=8 -> 9=1-8
9-2=4 -> 6-2=4
9-2=5 -> 3+2=5
9-2=6 -> 9-3=6
9-2=7 -> 9=2-7
9-2=7 -> 5+2=7
9-2=8 -> 8-2=6
9-3=0 -> 9-3=6
9-3=3 -> 6-3=3
9-3=4 -> 9-5=4
9-3=6 -> 9=3-6
9-3=6 -> 3+3=6
9-3=7 -> 9-2=7
9-3=8 -> 5+3=8
9-3=8 -> 9-9=0
9-3=9 -> 9-3=6
9-3=9 -> 8-3=5
9-4=2 -> 6-4=2
9-4=3 -> 9-4=5
9-4=5 -> 9=4-5
9-4=7 -> 3+4=7
9-4=9 -> 5+4=9
9-5=1 -> 6-5=1
9-5=4 -> 9=5-4
9-5=6 -> 9-3=6
9-5=8 -> 3+5=8
9-5=8 -> 9-9=0
9-5=9 -> 8-5=3
9-6=0 -> 6-6=0
9-6=0 -> 9-9=0
9-6=2 -> 9-6=3
9-6=3 -> 9=6-3
9-6=5 -> 9-6=3
9-6=7 -> 9-8=1
9-6=9 -> 9-0=9
9-6=9 -> 3+6=9
9-7=0 -> 9-1=8
9-7=2 -> 9=7-2
9-7=3 -> 9-7=2
9-7=6 -> 9-1=8
9-7=7 -> 8-1=7
9-7=7 -> 8-7=1
9-7=9 -> 9-1=8
9-8=1 -> 9=8-1
9-8=2 -> 8-6=2
9-8=3 -> 9-0=9
9-8=5 -> 9-0=9
9-8=8 -> 8-0=8
9-8=8 -> 8-8=0
9-8=9 -> 9+0=9
9-9=0 -> 9=9-0
9-9=3 -> 9-6=3
9-9=3 -> 8-5=3
9-9=5 -> 8-3=5
9-9=6 -> 9-9=0
9-9=7 -> 9-8=1
9-9=9 -> 9-0=9
9-9=9 -> 9-9=0

Mam nadzieję, że nie ma jakich drastycznych błędów w rozwiązaniu?

W ramach ćwiczeń polecam każdemu chętnemu rozwiązanie tego problemu w Waszym ulubionym języku programowania. Jeżeli twierdzicie, że da się prościej to zapewne macie racje, ale zawsze miło mi zobaczyć Wasz kod.

| Komentarze

Zdjęcia z kilku ostatnich lotów dronem

Drony pozwalają na spojrzenie na świat dookoła nas z innej perspektywy. W przeciągu ostatnich tygodni latałem kilkanaście razy i chciałem podzielić się z Wami zdjęciami, które udało mi się zrobić. Znad Wisły w okolicy mostu Świętokrzyskiego zrobiłem kila ujęć obydwu stron Warszawy:

dji_0010-hdrdji_0015-hdrdji_0069-hdr

Wyjątkowo ładny efekt udało mi się osiągnąć poprzez połączenie 5 ujęć (wszystkie HDR) w panoramę:

dji_0039-hdr-pano

Niewiele dalej na północ, przy moście Gdańskim zrobiłem kilka ujęć Warszawy w nocy:

dji_0005dji_0021dji_0049

Jeden dzień był wyjątkowo pochmurny, ale okazało się, że całkiem blisko świeci słońce:

dji_0195-hdrdji_0219-hdrdji_0226-hdr

Z lotu w chmury mam również film:

Jeżeli podobają się Wam moje zdjęcia i filmy to obserwujcie mnie na Instagramie (https://www.instagram.com/m.polkowski/) i subskrybujcie na YouTube (https://www.youtube.com/channel/UC0rcoI-FVLZjhIgTgfQpl9Q). Zapraszam!

| Komentarze

GMT – piękne mapy w kilku krokach

GMT (Generic Mapping Toolbox) to zestaw narzędzi do pracy z danymi geograficznymi (i nie tylko) opracowany na uniwersytecie na Hawajach (o mało co nie przyjęli mnie tam do pracy nad tym projektem… ehhh…), o którym mało kto pewnie słyszał. Podstawowym zastosowaniem GMT jest tworzenie różnego rodzaju map. Na stronie domowej projektu opisana jest instalacja (Windows, Mac OS X, Linux) oraz udostępniona jest obszerna dokumentacja, więc nie będę tych rzeczy tutaj powtarzał. Chciałem Was zachęcić kilkoma prostymi przykładami do zapoznania się z tym pakietem. Uwaga: to nie jest pełny tutorial, to tylko zachęta!

Przykład 1: Mapa bazowa Europy. GMT jest oczywiście pakietem bez interfejsu graficznego, więc wszystko odbywa się z wiersza poleceń:

gmt psbasemap -R-15/35/30/72 -Jl9/11/45/70/1:40000000 -B5g2:."Europa": -K -P -Y10 > Europa.ps
gmt pscoast -R -Di -Ia/0.03p,black -J -N1/0.5p,- -G#dde6d5 -S#a5bfdd -A0/0/4 -W0.5p,black -O -P >> Europa.ps
convert -geometry 2048x2048 -density 300 -trim Europa.ps Europa.png


Pierwszą komendą tworzymy mapę bazową o podanym wymiarze (-R), w projekcji azymutalnej równopowierzchniowej Lamberta (-Jl). GMT wynik (w postaci języka post script) wypisuje na standardowe wyjście. Drugą komendą rysujemy na mapie linię brzegową i granice państw. Opcja -D odpowiada za dokładność granic, -I za rysowanie rzek, -G i -S za kolory morza i lądów itd. (wszystko dokładnie opisane w dokumentacji). Trzecią komendą konwertujemy ps na png. Oto wynik:

europa

Przykład 2: Mapa bazowa Europy z naniesionymi punktami. Sama mapa bazowa jest ładna, ale mało użyteczna. Często chcemy na mapie pokazać jakieś lokalizacje. Na przykład lotniska. Ze strony http://openflights.org/data.html#airport pobrałem plik z listą lotnisk na świecie i przygotowałem mapę:

awk -F, '{print $8 "\t" $7}' airports.dat > airports.txt
gmt psbasemap -R-15/35/30/72 -Jl9/11/45/70/1:40000000 -B5g2:."Airports in Europe": -K -P -Y10 > Europa_air.ps
gmt pscoast -R -Di -J -N1/0.5p,- -G#dde6d5 -S#a5bfdd -A0/0/1 -W0.5p,black -O -P -K >> Europa_air.ps
gmt psxy airports.txt -Jl -O -R -Sc0.1c -W -G255/0/0 -K >> Europa_air.ps
convert -geometry 2048x2048 -density 300 -trim Europa_air.ps Europa_air.png


W pierwszej linii konwertują plik z danymi lotnisk na listę współrzędnych geograficznych. Drugą i trzecią linię już znamy. W czwartej linii rysuję na mapie czerwone kółeczka w każdej z zapisanych lokalizacji. Oto wynik:

europa_air

Przykład 3: Mapa topograficzna Polski. Jako przykład prezentacji rozkładu przestrzennego zmiennej wybrałem wysokość terenu nad poziomem morza. Miałem przygotowany wcześniej plik tekstowy zawierający informacje o wysokości n.p.m.:

screenshot_42

Teraz wystarczy tylko przeliczyć te dane na siatkę geograficzną i nanieść na mapę:

gmt surface DEM.xyz -R13.8/24.5/48.7/55 -I0.02/0.01 -GDEM.grd -T0.5 -C0.01
gmt psbasemap -R13.8/24.5/48.7/55 -Jl19.20/19.30/52/60/1:5000000 -B2g2:."Polska - Mapa Topo": -K -P -Y10 > PL.ps
gmt makecpt -CDEM_PL_wiki.cpt -T-10/2500/1 > tmp.cpt
gmt psscale -D3.1i/-.5i/7.0i/0.5ih -Ctmp.cpt -B250::/:m: -O -K -P >> PL.ps
gmt grdimage -Bg DEM.grd -R -J -O -Q -K -nn+a -Ctmp.cpt >> PL.ps
gmt pscoast -R  -Ia -Df -J -N1/0.5p,- -S#a5bfdd -A0/0/4 -W0.25p,black -K -O -P -L23.7/49.1/54/100+u -T14.2/54.5/1.5:: >> PL.ps
convert -geometry 2048x2048 -density 300 -trim PL.ps PL.png

W pierwszej linii przeliczamy dane w punktach na powierzchnię w siatce geograficznej i zapisujemy do pliku grd. W drugiej linii tworzymy mapę bazową dla Polski. Trzecia linia to przygotowanie skali kolorów na podstawie wzorca (DEM_PL_wiki.cpt, dużo wzorców tutaj). Czwarta linia to dodanie do mapy skali kolorów. Piąta linia to naniesienie na mapę naszej powierzchni w kolorach zgodnych ze skalą. Szósta linia to nałożenie na wszystko granic i zamalowanie wody na niebiesko. Oto efekt:

pl

GMT ma duże możliwości, niestety kosztem łatwości obsługi. Ważną zaletą jest tworzenie map w postaci wektorowej. GMT polecałbym szczególnie w dwóch przypadkach: po pierwsze gdy potrzebujemy bardzo wysokiej jakości np. do publikacji i po drugie w przypadku, gdy chcemy proces automatyzować np. generując indywidualne mapy dla użytkowników itp.

Gdybyście mieli jakieś pytania, piszcie – używam GMT od 5 lat więc może będę mógł pomóc.

| Komentarze

Sprytne omijanie limitu zapytań z adresu IP

Dawno dawno temu podjąłem się parsowania danych z pewnej strony www, na której obowiązywał limit ilości zapytań z jednego adresu IP. Działało to mniej więcej tak, że po kilku zapytaniach dany adres IP był blokowany na kilka minut. Ponieważ do pobrania było bardzo dużo rekordów, filtr ten skutecznie mnie blokował. Zacząłem szukać rozwiązania i po kilku iteracjach doszedłem do prostego i relatywnie taniego sposobu.

Na rynku istnieje wielu dostawców usługi zwanej Hosting SEO (do pozycjonowania) – od zwykłego hostingu różni się możliwością wykorzystania wielu adresów IP w ramach jednego konta. Adresy IP są stosunkowo drogie (jeśli potrzebujemy kilkadziesiąt), ale w hostingu SEO ta sama pula adresów dostępna jest dla wielu kont klientów. Za konto z pulą 32 adresów IP zapłacimy rzędu 20-30 zł za miesiąc. Po wykupieniu usługi i zalogowaniu się do panelu admina powinniśmy mieć dostęp do listy “naszych” adresów ip.

Zanim przejdziemy do rozwiązania, zróbmy prosty test. Spróbujmy z poziomu serwera sprawdzić adres IP. Zrobimy to za pomocą prostego skryptu php:

<?php
	$ch = curl_init('http://whatismyip.org/');
	curl_setopt($ch,CURLOPT_RETURNTRANSFER,TRUE);
	$myIp = curl_exec($ch);
	preg_match_all('/\d{1,3}\.\d{1,3}\.\d{1,3}\.\d{1,3}/',$myIp, $out);
	print_r($out);
?>

Skrypt uruchamiamy przez przeglądarkę (tani serwer raczej nie oferuje dostępu przez ssh). Wynik powinien wyglądać mniej więcej tak:

Array
(
    [0] => Array
        (
            [0] => 212.59.244.5
        )

)

Adres IP, który się wyświetli jest najprawdopodobniej przypisany domyślnie do naszego konta. Możemy przez panel administratora go zmienić i spróbować ponownie. Niestety takie rozwiązanie nie pomoże nam w szybkim zmienianiu adresów. Okazuje się, że wystarczy przy wywoływaniu CURL dodać informację, z którego interfejsu sieciowego serwera chcemy skorzystać. Na liście dostępnych dla mojego konta adresów IP były 32 pozycje:

ip

Użyjmy pierwszego adresu z listy:

<?php
	$ch = curl_init('http://whatismyip.org/');
	curl_setopt($ch,CURLOPT_RETURNTRANSFER,TRUE);
	curl_setopt($ch,CURLOPT_INTERFACE,"31.6.69.41");
	$myIp = curl_exec($ch);
	preg_match_all('/\d{1,3}\.\d{1,3}\.\d{1,3}\.\d{1,3}/',$myIp, $out);
	print_r($out);
?>

Otrzymujemy wynik:

Array
(
    [0] => Array
        (
            [0] => 31.6.69.41
        )

)

Działa! Za pomocą CURL możemy zdziałać cuda. Wystarczy tylko dopisać logikę, która losowo będzie przy każdym zapytaniu zmieniała adresy IP.

Jeśli ktoś jest dociekliwy to może odkryć jeszcze jedną rzecz. Wykupując najtańsze konto (z najmniejszą liczbą adresów IP) będziemy działali prawdopodobnie na tej samej maszynie co właściciele droższych pakietów. Możemy więc korzystać z całej puli adresów IP przypisanych do tego serwera, nawet jeśli nie są dla nas dostępne w panelu. U mnie działały np. takie adresy spoza listy: 31.6.69.42, 31.6.69.43, 31.6.69.45, 31.6.69.47 (więcej nie próbowałem). Użycie adresu nie przypisanego do danej maszyny spowoduje zwrócenie pustego stringa – nasze dane nie wiedziały jak wrócić.

Trzeba pamiętać jeszcze o małej małych zasobach serwera dla takiego konta hostingowego – jeśli zaczniemy generować za duże obciążenie prędzej czy później odezwie się admin z pretensjami.

Może komuś z Was się takie coś kiedyś przyda.

| Komentarze

Reprezentacja czasu w komputerze

Czas jest jednym z podstawowych “parametrów” świata, który na co dzień wielokrotnie zapisujemy (zapamiętujemy) i przeliczamy. Obecnie oznacza to, że jest przetwarzany przez komputery. Mogło by się wydawać, że sprawa jest prosta, ale od razu widać całą masę kłopotów. Pomysł na ten wpis wynika z bezsensownej zmiany czasu na zimowy w minioną niedzielę.

Po pierwsze strefy czasowe. Dopóki znajdujemy się w jednej strefie czasowej wszystko jest ok, potrafimy bez problemu obliczyć różnicę czasu pomiędzy dwiema datami. Do czasu: mamy przecież w Polsce dwie strefy czasowe – letnią i zimową. Jeżeli interesuje nas dokładne obliczenie czasu to musimy te informacje uwzględniać. Rozwiązaniem tego problemu może być zapamiętywanie i obliczanie dat i czasów w UTC (Universal Time Coordinated). Wiedzę na temat aktualnej w danym momencie i miejscu strefie czasowej najlepiej pobierać z systemu operacyjnego, który powinien to “ogarniać”. Nie warto tego rozgryzać samemu, bo ilość wyjątków jest powalająca. W tym miejscu polecam film The Problem with Time & Timezones.

Do tej pory pamiętam kartkówkę z geografii w liceum na której liczyliśmy czas podróży samolotem z San Francisco to Tokio (i z powrotem) na podstawie godzin przylotów i odlotów. Jedynym rozwiązaniem gwarantującym dobry wynik było przejście do czasu UTC.

Po drugie format zapisu danych. Komputer operuje na liczbach a ludzie na ciągach znaków. Istnieje wiele formatów zapisu daty i czasu w formie ciągu znaków, ale to jest tylko reprezentacja “dla ludzi”. W pamięci data i czas mogą być zapisane na wiele sposobów, z których trzy są moim zdaniem najważniejsze:

  • Możemy pamiętać oddzielnie jako liczby całkowite rok, miesiąc, dzień, godzinę, minutę, sekundę i milisekundę (albo nanosekundę itp.) + dodatkowo kod strefy czasowej (zapisany np. w postaci liczby minut w stosunku do UTC). Taki sposób jest dokładny i uniwersalny, ale mało wydajny i pamięciożerny. Ktoś jeszcze pamięta problem roku 2000?
  • Drugim sposobem jest zapisywanie epoki unix, czyli ilości sekund od 01-01-1970. Tu problemy są dwa: stosując int32 mamy bardzo ograniczony zakres dat i dokładność na poziomie sekundy.
  • Trzecim sposobem, z którym spotykam się na co dzień jest zapisywanie liczby zmiennoprzecinkowej określającej dni od roku 0. Z tej metody korzysta Matlab. Zdecydowaną zaletą jest możliwość opisania dowolnej daty w historii i przyszłości wszechświata oraz możliwość uwzględniania dokładności sub-sekundowej. Wada tego zapisu wynika ze zmiennej dokładności zapisu zmiennoprzecinkowego w zależności od wartości – pisałem o tym kilka tygodni temu.

Zdecydowanie trzeci sposób zapisu daty jest moim ulubionym – pozwala na bardzo szybkie operowanie na datach, bez utraty precyzji. Operacje na liczbach zmiennoprzecinkowych są całkowicie rutynowe, więc ciężko będzie znaleźć coś wydajniejszego.

Bardzo wiele zależy od języka programowania, w każdym jest inne podejście. Szczególnie w Pythonie widać ogrom komplikacji. Jakie są Wasze doświadczenia? Z którego formatu korzystacie najczęściej?

Jako bonus chciałbym Wam pokazać 3 animacje mojego autorstwa. Pierwsza jest dość przewidywalna bo pokazuje czas trwania dnia w zależności od dnia roku:

Dwie kolejne są już dużo ciekawsze, ponieważ pokazują czas wschodu i zachodu słońca w lokalnej strefie czasowej.

Widać wyraźnie zmianę czasu letniego na zimowy (w różnym czasie zależnie od regionu) i słabe dopasowanie strefy czasowej w niektórych regionach.

| Komentarze

Czy Python może być 876 razy szybszy?

Pythona zacząłem używać jakieś 5 lat temu. Na początku nie byłem przekonany, ale nie miałem wyboru. Z czasem się przekonałem i zrobiłem kilka projektów małych i dużych. Dodatkowo uczyłem podstaw programowania w pythonie przez 4 lata – do tej pory nie jestem przekonany czy jest to dobry język na początek (trochę za dużo wybacza), ale nie była to moja decyzja.

Wracając do tematu: python ma jedną poważną wadę – jest strasznie wolny. Dziś chcę wam pokazać na przykładzie, jak w prosty sposób można użyć funkcji napisanych w C/C++ do przeprowadzania obliczeń w programie w napisanym w pythonie. Python jest rewelacyjny jeśli chodzi o wczytywanie, zapisywanie i wizualizację danych, z kolei C/C++ jest bardzo wydajny. Połączenie wydaję się idealne i okazuje się całkiem proste.

Dla przykładu weźmy proste zadanie: generujemy listę N liczb losowych z zakresu od 0 do 1. Chcemy obliczyć wszystkie możliwe sumy 3 liczb z tej listy bez powtórzeń i policzyć ile z nich wynosi 1 z pewną ustaloną dokładnością. Oczywiście można ten problem optymalizować, ale dziś nie o tym. Chcemy porównać wydajność czystego pythona i pythona wspieranego przez C++.

W pierwszej kolejności napiszmy sobie funkcję w C++ (plik CALC.cpp):

#include "CALC.h"
#include 
#include 
#include 

void _CALC(double *DATA, int *RESULT, int N,  double treshold)
{
	printf("C++ start\n");

	RESULT[0] = 0;
	for (int a = 0; a < N; a++)
	{
		for (int b = a+1; b < N; b++)
		{
			for (int c = b+1; c < N; c++)
			{
				if (fabs(DATA[a] + DATA[b] + DATA[c] - 1) < treshold)
				{
					RESULT[0]++;
				}
			}
		}
	}

	printf("C++ finish\n");
}

W parze mamy plik nagłówkowy (CALC.h):

void _CALC(double *DATA, int *RESULT, int N, double treshold);

Powyższa funkcja robi dokładnie to co założyliśmy - zlicza ile kombinacji bez powtórzeń daje sumę 1 z dokładnością "treshold". Za chwilę powinno być jasne, czym są argumenty tej funkcji.

Naszym celem jest możliwość wywołania tej funkcji z poziomu pythona w taki sposób (program.py):

import numpy
import CALC
import time
start_time = time.time()

N = 1000
treshold = 1e-7

DATA = numpy.random.uniform(0, 1, size=N)
RESULT = numpy.zeros((1,), dtype=numpy.int)

start_time = time.time()
CALC.CALC(DATA, RESULT, N, treshold)
t1 = (time.time() - start_time)
print("C++:", RESULT[0], t1)


start_time = time.time()
count = 0;
for a in range(0,N):
	for b in range(a+1,N):
		for c in range(b+1,N):
			if abs(DATA[a]+DATA[b]+DATA[c]-1) < treshold:
				count =count + 1
t2 = (time.time() - start_time)
print("Python:", count, t2)

print("Speed:", t2/t1)

Dokładnie chodzi o linię

CALC.CALC(DATA, RESULT, N, treshold)

Potrzebujemy powiązania (interfejsu) pomiędzy pythonem i C++ (CALCmodule.cpp). Jest to kawałek kody, który z naszej funkcji C++ pozwala zrobić moduł dla pythona. Nie pisałem tego sam, złożyłem z kilku tutoriali. Najważniejsza jest funkcja CALC, gdzie definiujemy nasze zmienne i wywołujemy funkcję liczącą:

#include 
#define NPY_NO_DEPRECATED_API NPY_1_10_API_VERSION
#include 
#include "CALC.h"



struct module_state {
    PyObject *error;
};

#define GETSTATE(m) ((struct module_state*)PyModule_GetState(m))

static PyObject* CALC(PyObject* self, PyObject* args)
{
	PyArrayObject *data_obj;
	PyArrayObject *result_obj;
	double treshold;
	int N;

	if (!PyArg_ParseTuple(args, "OOid", &data_obj, &result_obj, &N, &treshold))
	{
		Py_INCREF(Py_None);
		return Py_None;
	}

	double *DATA = static_cast(PyArray_DATA(data_obj));
	int *RESULT = static_cast(PyArray_DATA(result_obj));

	_CALC(DATA, RESULT, N, treshold);

	Py_INCREF(Py_None);
	return Py_None;
}

static PyMethodDef CALCMethods[] = {
    {"CALC", CALC, METH_VARARGS, "..."},
    {NULL, NULL, 0, NULL}
};

static int CALC_traverse(PyObject *m, visitproc visit, void *arg) {
    Py_VISIT(GETSTATE(m)->error);
    return 0;
}

static int CALC_clear(PyObject *m) {
    Py_CLEAR(GETSTATE(m)->error);
    return 0;
}


static struct PyModuleDef moduledef = {
        PyModuleDef_HEAD_INIT,
        "CALC",
        NULL,
        sizeof(struct module_state),
        CALCMethods,
        NULL,
        CALC_traverse,
        CALC_clear,
        NULL
};

extern "C" PyObject * PyInit_CALC(void)
{
	PyObject *module = PyModule_Create(&moduledef);
	if (module == NULL)
        return NULL;
    struct module_state *st = GETSTATE(module);

    st->error = PyErr_NewException("CALC.Error", NULL, NULL);
    if (st->error == NULL) 
    {
        Py_DECREF(module);
        return NULL;
    }
	import_array();
	Py_INCREF(module);
    return module;

}

Na stronie https://docs.python.org/2/c-api/arg.html znajdują się dokładnie opisy typów danych które możemy przekazywać z pythona do C++ i z powrotem. W powyższym przykładzie "OOid" oznacza dwa obiekty, liczbę całkowitą i zmiennoprzecinkową (double). Obiektem może być lista, z której możemy czytać i do której pisać. Nazwa naszego modułu pojawia się w tym kodzie wielokrotnie - zmieniając trzeba pamiętać o wszystkich miejscach.

Ostatnim krokiem jest kompilacja modułu. Potrzebujemy skryptu do kompilacji (setup.py):

from distutils.core import setup, Extension
import numpy.distutils.misc_util
import os

os.environ["CC"] = "g++"
os.environ["CXX"] = "g++"

module1 = Extension('CALC', sources = ['CALCmodule.cpp', 'CALC.cpp'])

setup (name = 'CALC',
       version = '1.0',
       description = 'Package for calculating number sums',
       ext_modules = [module1],
	   include_dirs=numpy.distutils.misc_util.get_numpy_include_dirs())

Teraz możemy kompilować i uruchamiać:

screenshot_35

Funkcja w C++ potrzebowała 0.12 sekundy na policzenie, że jest 19 trójek, które z dokładnością 0.0000001 dają sumę 1. Python potrzebował na analogiczne liczenie 108.42 sekund: 876 razy dłużej! Oczywiście można by ten kod zoptymalizować, ale porównujemy dwa identyczne rozwiązania. Różnica w czasie jest kolosalna. W mojej pracy często spotykam się z obliczeniami numerycznymi, dlatego w pierwszej kolejności przygotowuję prototyp w pythonie i fragment obliczeniowy przepisuję na C++.

Nie chcę wnikać w teorię budowania modułów dla pythona bo się na tym nie znam. Sam najlepiej uczę się na przykładach, więc mam nadzieję, że ten Wam się przyda.

Dajcie znać, czy kiedyś używaliście takich rozwiązań, albo czy przydadzą Wam się w przyszłości.

| Komentarze (2)

Starsze wpisy »