猫史档案馆

import time
import matplotlib.pyplot as plt

__all__ = ['fract', 'poly']

class fract:
    def __init__(self,*arg):
        arg=list(arg)
        if type(arg[0])==list:
            arg=arg[0]
        elif len(arg)==1:
            arg.append(1)
        self.a=arg[0]
        self.b=arg[1]
        self.reduce_fract()
    
    def __str__(self):
        if self.b==1:
            return f'{int(self.a)}'
        return f'{int(self.a)}/{int(self.b)}'

    def check(obj):
        if type(obj)==int:
            return fract(obj,1)
        return obj

    def __add__(self,other):
        other=fract.check(other)
        if self.b==other.b:
            b=self.b
            a=self.a+other.a
        else:
            b=self.b*other.b
            a=self.a*other.b+other.a*self.b
        return fract(a,b)
    
    def __sub__(self,other):
        return self+(-other)
    
    def __mul__(self,other):
        other=fract.check(other)
        a=self.a*other.a
        b=self.b*other.b
        return fract(a,b)
    
    def __truediv__(self,other):
        other=fract.check(other)
        a=self.a*other.b
        b=self.b*other.a
        return fract(a,b)
    
    def __neg__(self):
        return fract(-self.a,self.b)
    
    __radd__=__add__
    __rmul__=__mul__
    __rsub__=lambda self,other:other+(-self)
    
    def factor(x, y):
        y1=y
        x,y=int(abs(x)),int(abs(y))
        if max(x,y)%min(x,y)==0:
            i = min(x,y)
        else:
            for i in range(min(x, y)//2+1,0,-1):
                if x % i == 0 and y % i == 0:
                    break
        if y1 < 0:
            return -i
        else:
            return i

    def reduce_fract(self):
        if self.a==0:
            self.b=1
        else:
            f=fract.factor(self.a,self.b)
            self.a/=f
            self.b/=f



class poly:
    def __init__(self,*arg,is_list=False,down=None): # *arg:降幂 list:升幂
        if down is None:
            down = not is_list
        if is_list:
            item = arg[0]
        else:
            item = list(arg)
        if down:
            item = item[::-1]
        for i in range(len(item)):
            v = item[i]
            if type(v) in [int,list]:
                item[i] = fract(v)
            elif type(v) != fract:
                raise ValueError('参数必须为 int,list 或 fract')
        self.item=item
        self.reduce_poly()

    def output(self,ret=False):
        self.reduce_poly()
        item=self.item
        txt='f(x)='
        for i in range(len(item))[::-1]:
            a,b=item[i].a,item[i].b
            if a % 1 == 0:
                a = int(a)
            if b % 1 == 0:
                b = int(b)
            if a == 0:
                continue
            if i<len(item)-1 and a>0:
                txt+='+'
            if a==1:
                if i==0 or b!=1:
                    txt+='1'
            elif a == -1:
                if i == 0 or b != 1:
                    txt += '-1'
                else:
                    txt += '-'
            else:
                txt+=f'{a}'
            if b!=1:
                txt += f'/{b}'
            if i > 0:
                txt+='x'
            if i>1:
                txt+=f'^{i}'
        if ret:
            return txt
        print(txt)
    
    def __str__(self):
        return self.output(True)

    def __call__(self,x):
        item=self.item
        val=0
        for i in range(len(item)):
            val += x**i*item[i].a/item[i].b
        if val % 1 == 0:
            val=int(val)
        return val

    def print_values(self, x1=None, x2=None):
        lag = hasattr(self, 'lagx')
        if lag:
            if x1 is None:
                x1 = min(self.lagx)
            if x2 is None:
                x2 = max(self.lagx)
        for x in range(x1, x2+1):
            v=self(x)
            if type(v)==float:
                v = round(self(x),4)
                if v % 1 == 0:
                    v = int(v)
            print(f'f({x})={v}')

    def draw(self,x1=None,x2=None,*args,point=2):
        if hasattr(self, 'lagx'):
            if x1 is None:
                x1=min(self.lagx)
            if x2 is None:
                x2=max(self.lagx)
        xlen=x2-x1
        xs = [x1+xlen*x/100 for x in range(101)]
        for ply in [self]+list(args):
            ys = [ply(x) for x in xs]
            plt.plot(xs, ys, "r")
            if point == 1 or not hasattr(ply, 'lagx'):
                xl = [x for x in range(x1, x2+1)]
                yl = [ply(x) for x in xl]
                plt.plot(xl, yl, '.k ')
            elif point == 2:
                xl = [x[0]/x[1] if type(x) == list else x for x in ply.lagx if x1 <= x <= x2]
                yl = [ply(x) for x in xl]
                plt.plot(xl, yl, '.k ')
            if min(ys) <= 0 <= max(ys):
                plt.plot([x1, x2], [0, 0], 'k')
            if x1<=0<=x2:
                plt.plot([0,0],[min(ys), max(ys)], 'k')
        plt.show()

    def reduce_poly(self):
        item=self.item[::-1]
        for i in range(len(item)):
            if item[0].a!=0:
                break
            del item[0]
        for i in range(len(item)):
            item[i].reduce_fract()
        self.item = item[::-1]

    def item(obj):
        if type(obj)==poly:
            return obj.item[:]
        return [obj][:]

    def __add__(self, other):
        item1=poly.item(self)
        item2=poly.item(other)
        itemM = max(item1,item2,key=len)
        item=itemM[:]
        for i in range(len(item)):
            try:
                item[i]=item1[i]+item2[i]
            except IndexError:
                break
        return poly(item,is_list=True)
    
    def __sub__(self, other):
        return self+(-other)

    def __mul__(self,other):
        item1=poly.item(self)
        item2=poly.item(other)
        item = [fract(0) for i in range(len(item1)+len(item2))]
        for i1 in range(len(item1)):
            for i2 in range(len(item2)):
                it1, it2 = item1[i1], item2[i2]
                item[i1+i2] += it1*it2
        return poly(item,is_list=True)
    
    def __truediv__(self,frac):
        frac=fract.check(frac)
        if frac.a == 0:
            raise ZeroDivisionError('division by zero')
        item = self.item[:]
        for i in range(len(item)):
            item[i]/=frac
        return poly(item, is_list=True)
    
    def __neg__(self):
        item=self.item[:]
        for i in range(len(item)):
            item[i]=-item[i]
        return poly(item,is_list=True)

    __radd__=__add__
    __rmul__=__mul__
    __rsub__=lambda self,other:other+(-self)

    def lag(xl,yl=None,timing=False):
        t=time.time()
        if yl is None:
            yl=xl
            xl=list(range(1,len(yl)+1))
        xl,yl=xl[:],yl[:]
        ply = poly(0)
        for i in range(len(xl)):
            item=poly(yl[i])
            for j in range(len(xl)):
                if j!=i:
                    xi=fract(xl[i])
                    xj=fract(xl[j])
                    z=fract(0)
                    item *= poly(1,z-xj)
                    item /= xi-xj
            ply += item
        ply.lagx,ply.lagy=xl,yl
        if timing:
            t=round(time.time()-t,5)
            print(f'计算用时{t}s')
        return ply
    
    def lag114514(yl):
        return poly.lag(yl+[114514])

    def fit(func,xl):
        return poly.lag(xl,[func(x) for x in x喵14]
ply=lag(yl)
print(ply)

举例:1，3，5，7的下一项是什么?

答:是114514.根据田所浩二定理，易得当f(x)=114505/24x^4-572525/12x^3+4007675/24x^2-2862601/12x+114504
f(1)=1
f(2)=3
f(3)=5
f(4)=7
f(5)=114514

（确信）

ply=poly.lag(yl)