<jupyter_start><jupyter_text>**Authors:** <NAME>, <NAME>   *[Faculty of Science](https://www.upjs.sk/en/faculty-of-science/?prefferedLang=EN), P. J. Šafárik University in Košice, Slovakia*  emails: [<EMAIL>](mailto:<EMAIL>)***  Generating quadruple-precision $\mathcal{GDD}$ pdf values (PARI/GP) Computational tools - CAS:   **PARI/GP**    $\mathcal{GDD}$ pdf  $-$ analytic expression$f(z)=\dfrac{\beta_{1}^{\alpha_{1}} \beta_{2}^{\alpha_{2}}}{\beta^{\alpha-1}}\begin{cases}{\dfrac{e^{z \beta_{2}}}{\Gamma\left(\alpha_{2}\right)} U\left(1-\alpha_{2}, 2-\alpha,-z \beta\right),} & {z<0} \\[12pt]\begin{array}{cc} \frac{\Gamma(\alpha-1)}{\beta^{a / 2-1} \Gamma\left(\alpha_{1}\right) \Gamma\left(\alpha_{2}\right)}, & \scriptstyle 1<\alpha \\ \infty,  & \scriptstyle 0<\alpha \leq 1 \end{array}, & z=0 \\[12pt] {\dfrac{e^{-z \beta_{1}}}{\Gamma\left(\alpha_{1}\right)} U\left(1-\alpha_{1}, 2-\alpha, z \beta\right),} & {z>0}\end{cases}$$\alpha=\alpha_{1}+\alpha_{2}, \quad \beta=\beta_{1}+\beta_{2}$  PARI/GP code<jupyter_code>%%gp 

\p33        
\\ defining U and G
U(u,v,t) = hyperu(u,v,t);
G(s) = gamma(s);

[a1, b1, a2, b2] = [0.5, 1, 8.5, 93];
a = a1+a2; b = b1+b2;
c = b1^a1*b2^a2/b^(a-1);
cp = c/G(a1); cm = c/G(a2);

\\ defining f(z)
fm(z) = cm*exp( z*b2)*U(1-a2,2-a,-b*z);
f0(z) = if(a>1, c*G(a-1)/(G(a1)*G(a2)), oo);
fp(z) = cp*exp(-z*b1)*U(1-a1,2-a, b*z);

f(z) = if(z>0, fp(z), /*if*/ z<0, fm(z), /*else*/f0(z));

\\ generating function values
fvalues(n, dh, hh) = {flist = vector(n); for(i= 1, n, flist[i] = f(dh+(hh-dh)/(n-1)*(i-1))); flist};<jupyter_output><jupyter_text>Sage procedures and functions<jupyter_code># python libraries
import numpy as np
from numpy import array as v
from io import StringIO
from time import time

import matplotlib.pyplot as plt
plt.style.use('seaborn-whitegrid')
# approximate formulas for precisions expressed in bits and decimal places
bits = lambda d:round((d+1)*ln(10)/ln(2)) 
dps = lambda b:round(b*ln(2)/ln(10) - 1)
bits(33), dps(113)
# writing PARI pdf data in a file via print commnand
def write_pari(filename, n, dh, hh):
    f = open(filename,"w")
    data = gp("fvalues"+str((n,dh,hh)))
    datastring = str(data)[1:-1].replace(' ', '').replace('+oo','inf')
    print(datastring, file =f)
    f.close()

# reading PARI pdf data with quadruple precision
def read_pari(filename, delim = ',', dps = 33):
    RRf = RealField(round((dps+1)*ln(10)/ln(2)) )
    f = open(filename, "r")
    datastring = f.read()
    f.close()
    return [RRf(item) for item in datastring.replace('\n','').split(delim)]    

# generating PARI pdf with quadruple precision
def paridata(n, dh, hh, dps = 33):
    RRf = RealField(round((dps+1)*ln(10)/ln(2)) )
    data = gp("fvalues"+str((n,dh,hh)))
    datastring = str(data)[1:-1].replace(' ', '').replace('+oo','inf')
    return [RRf(item) for item in datastring.split(',')]

# writing x values with quadruple precision
def write_xdata(filename, n, dh, hh, dps = 33):
    f = open(filename,"w")
    RRf = RealField(round((dps+1)*ln(10)/ln(2)) )
    data = [RRf(QQ(dh)+(QQ(hh)-QQ(dh))/(n-1)*(i-1)) for i in [1..n]]
    for item in data:
        print(item, file=f)
    f.close()<jupyter_output><empty_output><jupyter_text>Testing and examples<jupyter_code># writing x values (quadruple precision) 
for n in [1..5]:
    tic = time()   # start time
    
    filename = 'xvalues'+str(10^n)+'.txt'
    write_xdata(filename, 10^n,-3,4)
    
    toc = time()-tic; # end time
    print('10^'+str(n)+' x values: runtime =',toc,'s')
# pre-calculated x values (quadruple precision) 
N = 5
dx = {str(10**(n+1)):np.loadtxt('xvalues'+str(10**(n+1))+'.txt', dtype=np.longdouble) for n in range(N)}
# writing 10^n pdf values in files
for n in [1..5]:
    tic = time()   # start time
    
    filename = 'Pari_Sage_pdf'+str(10^n)+'.txt'
    write_pari(filename, 10^n,-3,4)
    
    toc = time()-tic; # end time
    print('10^'+str(n)+' values: runtime =',toc,'s')
# plots
for n in [1..5]:
    filename = 'Pari_Sage_pdf'+str(10^n)+'.txt'
    data = read_pari(filename)
    g = line(zip(dx[str(10^n)],data))
    g.show(frame=True, title = '$\\mathcal{GDD}(0.5, 1, 8.5, 93)$ pdf $-$'+str(10^n)+' points')    
# generating 10 pdf values with PARIGP
n = 10
dh, hh = -3, 4 

tic = time()   # start time
data = gp("fvalues"+str((n,dh,hh)))
toc = time()-tic; # end time

toc
str(data), str(data)[1:-1]
# data conversion to one string
datastring = str(data)[1:-1].replace(' ', '').replace('+oo','inf')

# conversion to quadruple precision
RRf = RealField(113)
dparis = [RRf(item) for item in datastring.split(',')]
dparis
# read string via numpy loadtxt
dtxt = StringIO(datastring)
df = np.loadtxt(dtxt, delimiter=',', dtype = np.longdouble)
df
# save as csv
np.savetxt('test_Pari_Sage_pdf_10.txt', df)
# read file via numpy loadtxt
data = np.loadtxt('test_Pari_Sage_CDF_10.txt', delimiter=',')
data<jupyter_output><empty_output>
{'fraction_non_alphanumeric': 0.22199467540446446, 'fraction_numerical': 0.03604341593282818, 'mean_word_length': 4.240343347639485, 'pattern_counts': {'":': 0, '<': 20, '<?xml version=': 0, '>': 19, 'https://': 1, 'lorem ipsum': 0, 'www.': 1, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 0, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'chain_length': '3', 'dir': 'jupyter-structured-clean-dedup', 'id': '169008', 'license': 'MIT', 'n_tokens_mistral': 2123, 'n_tokens_neox': 1994, 'n_words': 436, 'path': 'Jupyter/3.2 Quadruple-precision pdf - PARI.ipynb', 'repo_name': 'fdslrm/GDD', 'size': '4928'}