Προγραμματιστική πρόκληση Νο 5: Υπολογισμός της τετραγωνικής ρίζας

Ελπίζω να μην είναι κανένας που να μην ξέρει τι είναι η τετραγωνική ρίζα ενός αριθμού. Και σίγουρα ξέρει πως να την βρει στην αγαπημένη του γλώσσα προγραμματισμού. Αλλά μπορεί αλήθεια να βρεθεί χωρίς την αυτή κάποιας έτοιμης συνάρτησης ή τελεστή (ή την ύψωση σε δύναμη που είναι το ίδιο πράγμα);

Με άλλα λόγια πως θα υπολογίσεις μια τετραγωνική ρίζα, αν το μόνο που έχεις είναι μολύβι και χαρτί; Έστω ένα απλό κομπιουτεράκι, μιας και με την σκέψη μόνο πιθανά να ανατρίχιασες :smiley:

Για σκεφτείτε λίγο πως θα την υπολογίσετε; Δείχνω δυο τρόπους αλλά περιμένω την ιδέα σας πρώτα πριν τα ανοίξετε

Μέθοδος 1

Ας ξεκινήσουμε με μια προσεγγιστική τιμή. Την πολλαπλασιάζουμε με τον εαυτό της, και ανάλογα την αυξάνουμε ή την μειώνουμε. Αφήνω σε εσάς τις λεπτομέρειες.

Μέθοδος 2

Θα αφήσω τις βαρετές λεπτομέρειες και θα μπω στο ζουμί. Πάρτε ένα κομπιουτεράκι και βάλτε ένα τυχαίο αριθμό. Βρείτε το συνημίτονο αυτού του αριθμού (cos) και μετά το συνημίτονο αυτού κλπ μέχρι να μην αλλάζει άλλο ο αριθμός. Μπορείτε να γράψετε ένα πρόγραμμα και να βρείτε ποιος αριθμός είναι αυτός;

Τον αριθμό αυτό θα τον πούμε το σταθερό σημείο της συνάρτησης \cos(x). Εύκολα(*) αποδεικνύετε πως ο αριθμός \sqrt{\alpha} είναι το σταθερό σημείο της συνάρτησης f(x) = \frac{a}{x}.

Υπάρχει και η Βαβυλωνιακή μεθοδος ή μέθοδος του Ηρωνα. Με μια λογική παρόμοια με την πρώτη μέθοδο αν το x_n είναι μια προσέγγιση, τότε ο αριθμός x_{x+1} = \frac{1}{2}\left(x_n + \frac{\alpha}{x_n}\right) είναι μια καλύτερη προσέγγιση.

Ας το κάνουμε σε ταληράκια, και ας πούμε πως η πρώτη προσέγγιση για την τετραγωνική ρίζα είναι ο αριθμός 1. Με βάση τον τύπο θα έχουμε όλο και καλύτερες προσεγγίσεις της τετραγωνικής ρίζας. Μέχρι που θα είμαστε ικανοποιημένοι και θα σταματήσουμε. Όμορφο;

(*) Kωδική φράση των μαθηματικών όταν υπονοούν πως μετά από 3 ώρες πράξεων θα βρεις το αποτέλεσμα, οπότε καλύτερα να το βουλώσει. Αν δεν το βουλώσει θα ανακαλύψει το μυστικό πίσω την ίδια την έννοια του υπολογισμού, αλλά αυτό είναι μια μεγάλη ιστορία. Αν κάποιος θέλει λεπτομέρειες ας ξεκινήσει από εδώ.

Άλλες μέθοδοι και αλγόριθμοι μπορούν να βρεθούν εδώ.

Φτιάξτε μια συνάρτησή που να τις δίνουμε δυο εισόδους, ένα θετικό αριθμό και το πλήθος των δεκαδικών ψηφίων που θέλουμε και να υπολογίζει την τετραγωνική ρίζα του αριθμού με την ζητούμενη ακρίβεια. Πόσο γρήγορος είναι ο αλγόριθμος σας;

1 Like

Απαντώ έτσι πρόχειρα, στο πόδι που λένε, με τον καφέ. Δεν ξέρω αν είναι επιθυμητή η απάντηση. Είναι όμως χωρίς βιβλιοθήκη. :slight_smile:


num = input("Number: ")
result = int(num)**0.5
print(f'Η τετραγωνική ρίζα του {num} είναι {result}')

Η ύψωση σε κλασματική δύναμη που κάνεις είναι χρήση μιας έτοιμης βιβλιοθήκης που την χρησιμοποιείς έμμεσα :smiley:

Update Επίσης δεν λύνει το πρόγραμμα σου το ζητούμενο. Ας το δούμε λίγο ποιο προσεκτικά, δεν ζητάει τον υπολογισμό απλά μιας τετραγωνικής ρίζας.

Εσύ απλά υπολογίζεις έναν αριθμό, χωρίς κάποια εγγύηση ή πρόβλεψη για την δεκαδική ακρίβεια του υπολογισμού. Ναι υποθέτω θα είναι καλή ακρίβεια σε κάποιο Pentium, πολύ μεγαλύτερη από (πόσα αλήθεια;) τα 3-4 που αρκούν στους υπολογισμούς συνήθως, αλλά δεν είναι αυτό το ερώτημα.

Την math την έβαλα για έναν έλεγχο. Δεν την χρησιμοποιώ κάπου. Θεωρείται χρήση βιβλιοθήκης η ύψωση σε δύναμη;; Δεν το ήξερα. Τότε οποιαδήποτε πράξη θεωρείται χρήση βιβλιοθήκης, άρα όλο είναι άτοπο :joy:

1 Like

Τεχνικά η ύψωση σε δύναμη και η ρίζες είναι η χρήση μιας συνάρτησης βιβλιοθήκης. Συνήθως η βιβλιοθήκη αυτή είναι μια Intrinsic function και θα μεταφραστεί με ειδικό τρόπο σε μια εντολή της FPU του επεξεργαστή. Αλλά δεν έχουν όλοι οι επεξεργαστές FPU ή ακόμα και αν έχουν μπορεί να έχουν δυνατότητες μόνο πρόσθεση και πολλαπλασιασμό (μπορεί μαζί σε μία πράξη μόνο). Αυτό αρκεί για να υπολογίσεις ρίζες τριγωνομετρικές συναρτήσεις και άλλα. Αν δεν έχεις την εντολή θα πρέπει να γράψει κάποιος αυτόν των κώδικα.

O πυρήνας του Linux ίσως ακόμα να έχει ένα module που εξομοιώνει την FPU για επεξεργαστές που δεν έχουν (Για 486SX έχει αφαιρεθεί, για τους ARM υπήρχε μια περυσινή πρόταση να αφαιρεθεί, αλλά δεν ξέρω αν υλοποιήθηκε, αλλά ίσως υπάρχει ακόμα για άλλες αρχιτεκτονικές).

Αλλά πως η FPU θα κάνει τον υπολογισμό; Πάλι κάποιος άλλος θα πρέπει να γράψει αυτόν το κώδικά σε ειδικές γλώσσες όπως η Verilog και η VHDL.

Γιατί να κάνεις κάποιος άλλος αυτή την δουλεία; Μην ξεχνάμε πως η σειρά αυτή είναι μια σειρά από προκλήσεις απέναντι στον εαυτό σου :innocent: Σαν τον Ευκλείδη ας υποθέσουμε πως αντι για κανόνα και διαβήτη διαθέτεις τις βασικές πράξεις για ακεραίους και αριθμούς κινητής υποδιαστολής.

Θα κλείσω με ένα κώδικα που έχει κάνει μάτια και μυαλά να δακρύσουν σε προγραμματιστικές συνεντεύξεις. Τι λέτε να κάνει;

float Q_rsqrt( float number )
{
	long i;
	float x2, y;
	const float threehalfs = 1.5F;

	x2 = number * 0.5F;
	y  = number;
	i  = * ( long * ) &y;                       // evil floating point bit level hacking
	i  = 0x5f3759df - ( i >> 1 );               // what the fuck? 
	y  = * ( float * ) &i;
	y  = y * ( threehalfs - ( x2 * y * y ) );   // 1st iteration
//	y  = y * ( threehalfs - ( x2 * y * y ) );   // 2nd iteration, this can be removed

	return y;
}

Όχι δεν ήταν το θέμα της ερώτησης να πεις αυτό. Ήταν να δούνε πόσο θα κομπλάρεις :stuck_out_tongue: Είναι από τον κώδικα του Quake. Κάποιες φορές πρέπει να κατέβουμε χαμηλά. Για περισσότερα εδώ.

Αλλά εδώ είπαμε, βάζουμε προκλήσεις στους εαυτούς μας, είτε κάνουμε προγραμματισμό δεκαετίες είτε ξέρουμε ελάχιστο, για να ξεσκουρίαζουμε και να μαθαίνουμε :smiley_cat: Και γιατί να υλοποιήσεις μόνο μια λύση και όχι πολλές και να τις συγκρίνεις :grinning:

Χωρίς ύψωση σε δύναμη, με ακρίβεια 2 δεκαδικών

number = float(input("Δώσε αριθμό για υπολογισμό τετραγωνικής ρίζα : "))
riza_under = 0
step = 0.001
delta = 0.0001
while riza_under * riza_under < number and abs(riza_under * riza_under - number) > delta:
    riza_under = riza_under + step

riza_over = riza_under + step
if abs(riza_over * riza_over - number) < abs(riza_under * riza_under - number):
    print("{:.2f}".format(riza_over))
else:
    print("{:.2f}".format(riza_under))

Αν θες να μου πεις για τη συνάρτηση abs() προφανώς είναι το μόνο εύκολο να κάνω μια συνάρτηση που να κάνει ακριβώς το ίδιο.

2 Likes

Ας ρίξουμε μια ματιά και σε αυτόν τον κώδικα για τον υπολογισμό της τετραγωνικής ρίζας (δεν τον σκέφτηκα μόνος μου).Μου φαίνεται αρκετά ενδιαφέρον.

def sroot(n):
  sqrt = n/2
  temp = 0
  while sqrt != temp:
    temp = sqrt
    sqrt = ( n/temp + temp) / 2
  return sqrt
2 Likes

Ο κώδικάς είναι κατά βάση σωστός. Και στην εκφώνηση δίνω και ένα όνομα σε αυτόν τον αλγόριθμο. Αλλά για να δούμε λίγο καλύτερα το ζητούμενο

Η λύση σου παίρνει μόνο μια είσοδο. Επίσης, τι θα συμβεί αν δώσεις αρνητικό αριθμό;

Εγώ θα αυτό θα το έλεγα σοβαρό λάθος, όχι απλά ότι αγνοεί την εκφώνηση. Προσωπικά αν έφτιαχνα δική μου γλώσσα, η χρήση == ή != με κινητής υποδιαστολής θα πέταγε το λιγότερο warning.

Αυτό εγώ δεν το θεωρώ λάθος γιατί μπορείς συνήθως να ασχοληθείς με τα προβληματικά δεδομένα έξω από τον αλγόριθμο, ακόμα περισσότερο, μάλιστα το θεωρώ καλύτερη τακτική.

Άντε γκρινιάρηδες το έκανα λίγο καλύτερο

number = float(input("Δώσε αριθμό για υπολογισμό τετραγωνικής ρίζα : "))
precision = input("Δώσε ακρίβεια δεκαδικού ψηφίου : ")
riza_under = 0
step = 0.001
while riza_under * riza_under < number:
    riza_under = riza_under + step

riza_over = riza_under + step
if abs(riza_over * riza_over - number) < abs(riza_under * riza_under - number):
    print("{:.{}f}".format(riza_over,precision))
else:
    print("{:.{}f}".format(riza_under,precision))
2 Likes

Εγώ έλεγα για τον Χαράλαμπο

1 Like

Όπως προαναφέρθηκε, με τη μέθοδο του Ήρωνα, που είναι απλουστευμένη μορφή της μεθόδου Newton–Raphson για εύρεση της ρίζας, στην αριθμητική ανάλυση, υλοποίηση σε FORTRAN με δυνατότητα εκτύπωσης μέχρι 10 δεκαδικών ψηφίων και υπολογισμού μέχρι τον ακέραιο αριθμό 2147483647(που είναι ο μεγαλύτερος απλός ακέραιος στην συγκεκριμένη γλώσσα), μετά χρειάζεται error handling :

PROGRAM newtonsqrt
IMPLICIT NONE
INTEGER :: xx,n
PRINT*,"Dwse enan arithmo : "
READ*,xx
IF (xx<0)THEN
	PRINT*,"Dwsate arnitiko arithmo!"
	STOP
END IF
PRINT*,"Dwse arithmo dekadikwn (mexri 10 psifia): "
READ*,n
IF ((n>10).OR. (n<0))THEN
	PRINT*,"Den yparxei tetoia dinatotita!"
	STOP
END IF
CALL newtons(xx,n)
END PROGRAM newtonsqrt

SUBROUTINE newtons(xx,n)
IMPLICIT NONE
INTEGER :: xx,n,i
REAL(8) :: x
CHARACTER (LEN=21) :: stringx,tempstring
n=n+10
x=1
DO i=1,100
	x=0.5*(x+(xx/x))
	WRITE(stringx,'(F20.10)')x
	IF(stringx(1:n)==tempstring(1:n))THEN
		PRINT*,"Meta apo ",i," epanalipseis to apotelesma einai : ",stringx(1:n)
		RETURN
	END IF
	tempstring=stringx
END DO
RETURN
END SUBROUTINE newtons
2 Likes

Να προσθέσω και την μέθοδο του Bakhshali που έχει τετραπλάσια ταχύτητα στην προσέγγιση της λύσης :

SUBROUTINE bakh(xx,n)
IMPLICIT NONE
INTEGER :: i,n,xx
REAL (8):: x0,x,a1,b1
CHARACTER (LEN=21) :: stringx,tempstring
n=n+10
i=1
x0=1
DO i=1,20
	a1=(xx-x0*x0)/(2*x0)
	b1=x0+a1
	x=b1-(a1*a1)/(2*b1)
	WRITE(stringx,'(F20.10)')x
	IF(stringx(1:n)==tempstring(1:n))THEN
		PRINT*,"Meta apo ",i," epanalipseis to apotelesma einai : ",stringx(1:n)
		RETURN
	END IF
	PRINT*,"Meta apo ",i," epanalipseis to apotelesma einai : ",stringx(1:n)
	x0=x
	tempstring=stringx
END DO
END SUBROUTINE
1 Like