from __future__ import division
import math
def dypar(p,a1,m1,m2,e=1e-10): #p - perioda v letech, a1 - hlavni poloosa v uhl. vterinach, m1,2 - poz. hv. velikosti, e - presnost vypoctu
    n=0 #counter (zbytecnost - jen pro zajimavost a kontrolu poctu iteraci)
    mpc=3.082e16 #prepocet metru na pc
    M0=2e30 #solar mass
    M=[M0,M0] #vstupni hmotnost pro prvni odhad hl. poloosy
    m=[m1,m2] #zadane pozorovane hv. vel.
    Ps=p*365.25*86400 #prepocet periody na sekundy z let
    G=6.67e-11 #gravitacni konstanta
    a_ang=a1/3600/180*math.pi #prepocet uhlov0 paralaxy z vterin na radiany
    M_test=0 #hmotnost se kterou cyklus porovnava vysledek na ukonceni vypocty (nutna pouze pro prvni "prujezd" pak se bere predchozi)
    while abs(M_test-M[0])>e:
        M_test=M[0] # prirazeni porovn. hmotnosti pred vypoctem v nasl. cyklu
        a=(G*(M[0]+M[1])*Ps**2/4/math.pi**2)**(1/3) # vypocet hl. poloosy
        r=a/a_ang/mpc #vypocet vzdalenosti
        M_abs=[m[0]-5*math.log10(r)+5, m[1]-5*math.log10(r)+5] #vypocet abs. hv velikosti
        L=[10**((4.75-M_abs[0])/2.5),10**((4.75-M_abs[1])/2.5)] # vypocet L
        for i in range(2):  #prirazeni zavislosti podle M a spocteni nove M
            if M[i]<=0.43*M0:
                M[i]=(L[i]/0.23)**(1/2.3)*M0
            elif M[i]<=2*M0:
                M[i]=(L[i])**(1/3.5)*M0
            elif M[i]<=20*M0:
                M[i]=(L[i]/1.5)**(1/4)*M0
            else:
                M[i]=L[i]/3200*M0
        n+=1
    print ('Hmotnosti obou slozek v hmotnostech slunce: {0:.3f}, {1:.3f}.'.format(M[0]/M0,M[1]/M0))
    print ('Vzdalenost systemu: {0:.2f} pc'.format(r))
    print ('Dynamicka paralaxa: {0:.2f} mas'.format(1/r*1000))
    print ('Pocet iteraci ve vypoctu: {0}'.format(n))

dypar(87.85,4.551,3.93,5.29)