text:
from __future__ import print_function
import gc, os, sys, inspect
from copy import deepcopy
from warnings import warn
from time import time
import numpy as np
import scipy as sp
from numpy.linalg import norm
from Florence.Utils import insensitive
from Florence.Tensor import trace, Voigt

from Florence.PostProcessing import *
from Florence.Solver import LinearSolver, FEMSolver

__all__ = ["CoupleStressSolver"]

class CoupleStressSolver(FEMSolver):

    def __init__(self,**kwargs):
        if 'static_condensation' in kwargs.keys():
            self.static_condensation = kwargs['static_condensation']
        else:
            self.static_condensation = True

        if 'assembly_print_counter' in kwargs.keys():
            self.assembly_print_counter = kwargs['assembly_print_counter']
            del kwargs['assembly_print_counter']
        else:
            self.assembly_print_counter = 500

        if 'lump_rhs' in kwargs.keys():
            self.lump_rhs = kwargs['lump_rhs']
            del kwargs['lump_rhs']
        else:
            self.lump_rhs = False

        super(CoupleStressSolver, self).__init__(**kwargs)
        self.gamma   = 0.5
        self.beta    = 0.25

    def GetBoundaryInfo(self, mesh, formulation, boundary_condition):

        all_dofs = np.arange(mesh.points.shape[0]*formulation.nvar)
        if formulation.fields == "electro_mechanics" or formulation.fields == "flexoelectric":
            self.electric_dofs = all_dofs[formulation.nvar-1::formulation.nvar]
            self.mechanical_dofs = np.setdiff1d(all_dofs,self.electric_dofs)

            # GET BOUNDARY CONDITON FOR THE REDUCED MECHANICAL SYSTEM
            self.columns_in_mech = np.intersect1d(boundary_condition.columns_in,self.mechanical_dofs)
            self.columns_in_mech_idx = np.in1d(self.mechanical_dofs,boundary_condition.columns_in)

            # GET BOUNDARY CONDITON FOR THE REDUCED ELECTROSTATIC SYSTEM
            self.columns_in_electric = np.intersect1d(boundary_condition.columns_in,self.electric_dofs)
            self.columns_in_electric_idx = np.in1d(self.electric_dofs,boundary_condition.columns_in)

            # GET FREE MECHANICAL DOFs
            self.columns_out_mech = np.intersect1d(boundary_condition.columns_out,self.mechanical_dofs)
            self.columns_out_mech_idx = np.in1d(self.mechanical_dofs,boundary_condition.columns_out)

            # GET FREE ELECTROSTATIC DOFs
            self.columns_out_electric = np.intersect1d(boundary_condition.columns_out,self.electric_dofs)
            self.columns_out_electric_idx = np.in1d(self.electric_dofs,boundary_condition.columns_out)

            self.applied_dirichlet_mech = boundary_condition.applied_dirichlet[np.in1d(boundary_condition.columns_out,self.columns_out_mech)]
            self.applied_dirichlet_electric = boundary_condition.applied_dirichlet[np.in1d(boundary_condition.columns_out,self.columns_out_electric)]

            # MAPPED QUANTITIES
            out_idx = np.in1d(all_dofs,boundary_condition.columns_out)
            idx_electric = all_dofs[formulation.nvar-1::formulation.nvar]
            idx_mech = np.setdiff1d(all_dofs,idx_electric)

            self.all_electric_dofs = np.arange(mesh.points.shape[0])
            self.electric_out = self.all_electric_dofs[out_idx[idx_electric]]
            self.electric_in = np.setdiff1d(self.all_electric_dofs,self.electric_out)

            self.all_mech_dofs = np.arange(mesh.points.shape[0]*formulation.ndim)
            self.mech_out = self.all_mech_dofs[out_idx[idx_mech]]
            self.mech_in = np.setdiff1d(self.all_mech_dofs,self.mech_out)

            # LOCAL NUMBERING OF FIELDS
            self.all_local_dofs = np.arange(mesh.elements.shape[1]*formulation.nvar)
            self.all_local_electric_dofs = self.all_local_dofs[formulation.nvar-1::formulation.nvar]
            self.all_local_mech_dofs = np.setdiff1d(self.all_local_dofs,self.all_local_electric_dofs)

            # INTERMIX RAVELLED LOCAL INDICES
            matrix_shape = (mesh.elements.shape[1]*formulation.nvar,mesh.elements.shape[1]*formulation.nvar)
            self.idx_uu = np.ravel_multi_index(np.meshgrid(self.all_local_mech_dofs,self.all_local_mech_dofs), matrix_shape)
            self.idx_up = np.ravel_multi_index(np.meshgrid(self.all_local_mech_dofs,self.all_local_electric_dofs), matrix_shape)
            self.idx_pu = np.ravel_multi_index(np.meshgrid(self.all_local_electric_dofs,self.all_local_mech_dofs), matrix_shape)
            self.idx_pp = np.ravel_multi_index(np.meshgrid(self.all_local_electric_dofs,self.all_local_electric_dofs), matrix_shape)

        elif formulation.fields == "mechanics" or formulation.fields == "couple_stress":
            self.electric_dofs = []
            self.mechanical_dofs = all_dofs
            self.columns_out_mech = boundary_condition.columns_out

            self.mech_in = boundary_condition.columns_in
            self.mech_out = boundary_condition.columns_out

            self.applied_dirichlet_mech = boundary_condition.applied_dirichlet

    def __checkdata__(self, material, boundary_condition,
        formulation, mesh, function_spaces, solver, contact_formulation=None):

        ## IMPORTANT
        self.requires_geometry_update = False
        self.analysis_nature = "linear"

        function_spaces, solver = super(CoupleStressSolver,self).__checkdata__(material, boundary_condition,
            formulation, mesh, function_spaces, solver, contact_formulation=contact_formulation)

        # BUILD THE TANGENT OPERATORS BEFORE HAND
        if material.mtype == "CoupleStressModel":
            I = np.eye(material.ndim,material.ndim)

            material.elasticity_tensor = Voigt(material.lamb*np.einsum('ij,kl',I,I)+material.mu*(np.einsum('ik,jl',I,I)+np.einsum('il,jk',I,I)),1)
            material.gradient_elasticity_tensor = 2.*material.eta*I
            material.coupling_tensor0 = np.zeros((material.elasticity_tensor.shape[0],
                material.gradient_elasticity_tensor.shape[0]))

            # print material.elasticity_tensor
            ngauss = function_spaces[0].AllGauss.shape[0]
            d0 = material.elasticity_tensor.shape[0]
            d1 = material.elasticity_tensor.shape[1]
            d2 = material.gradient_elasticity_tensor.shape[0]
            d3 = material.gradient_elasticity_tensor.shape[1]
            material.elasticity_tensors = np.tile(material.elasticity_tensor.ravel(),ngauss).reshape(ngauss,d0,d1)
            material.gradient_elasticity_tensors = np.tile(material.gradient_elasticity_tensor.ravel(),ngauss).reshape(ngauss,d2,d3)

        elif material.mtype == "IsotropicLinearFlexoelectricModel":
            I = np.eye(material.ndim,material.ndim)

            material.elasticity_tensor = Voigt(material.lamb*np.einsum('ij,kl',I,I)+material.mu*(np.einsum('ik,jl',I,I)+np.einsum('il,jk',I,I)),1)
            material.gradient_elasticity_tensor = 2.*material.eta*I
            material.coupling_tensor0 = np.zeros((material.elasticity_tensor.shape[0],
                material.gradient_elasticity_tensor.shape[0]))
            material.piezoelectric_tensor = material.P
            material.flexoelectric_tensor = material.f
            material.dielectric_tensor = -material.eps*np.eye(material.ndim,material.ndim)

            # print material.elasticity_tensor
            ngauss = function_spaces[0].AllGauss.shape[0]
            d0 = material.elasticity_tensor.shape[0]
            d1 = material.elasticity_tensor.shape[1]
            d2 = material.gradient_elasticity_tensor.shape[0]
            d3 = material.gradient_elasticity_tensor.shape[1]
            material.elasticity_tensors = np.tile(material.elasticity_tensor.ravel(),ngauss).reshape(ngauss,d0,d1)
            material.gradient_elasticity_tensors = np.tile(material.gradient_elasticity_tensor.ravel(),ngauss).reshape(ngauss,d2,d3)
            material.piezoelectric_tensors = np.tile(material.P.ravel(),ngauss).reshape(ngauss,material.P.shape[0],material.P.shape[1])
            material.flexoelectric_tensors = np.tile(material.f.ravel(),ngauss).reshape(ngauss,material.f.shape[0],material.f.shape[1])
            material.dielectric_tensors = np.tile(material.dielectric_tensor.ravel(),ngauss).reshape(ngauss,material.ndim,material.ndim)

            factor = -1.
            material.H_Voigt = np.zeros((ngauss,material.H_VoigtSize,material.H_VoigtSize))
            for i in range(ngauss):
                H1 = np.concatenate((material.elasticity_tensor,factor*material.piezoelectric_tensor),axis=1)
                H2 = np.concatenate((factor*material.piezoelectric_tensor.T,material.dielectric_tensor),axis=1)
                H_Voigt = np.concatenate((H1,H2),axis=0)
                material.H_Voigt[i,:,:] = H_Voigt

        return function_spaces, solver

    def Solve(self, formulation=None, mesh=None,
        material=None, boundary_condition=None,
        function_spaces=None, solver=None,
        contact_formulation=None):
        """Main solution routine for FEMSolver """

        # CHECK DATA CONSISTENCY
        mesh = formulation.meshes[0]
        #---------------------------------------------------------------------------#
        function_spaces, solver = self.__checkdata__(material, boundary_condition,
            formulation, mesh, function_spaces, solver, contact_formulation=contact_formulation)
        #---------------------------------------------------------------------------#
        caller = inspect.getouterframes(inspect.currentframe(), 2)[1][3]
        if caller != "Solve":
            self.PrintPreAnalysisInfo(mesh, formulation)
        #---------------------------------------------------------------------------#

        # INITIATE DATA FOR NON-LINEAR ANALYSIS
        NodalForces, Residual = np.zeros((mesh.points.shape[0]*formulation.nvar,1),dtype=np.float64), \
            np.zeros((mesh.points.shape[0]*formulation.nvar,1),dtype=np.float64)
        # SET NON-LINEAR PARAMETERS
        self.NRConvergence = { 'Increment_'+str(Increment) : [] for Increment in range(self.number_of_load_increments) }

        # ALLOCATE FOR SOLUTION FIELDS
        TotalDisp = np.zeros((mesh.points.shape[0],formulation.nvar,self.number_of_load_increments),dtype=np.float64)
        TotalW = np.zeros((formulation.meshes[1].points.shape[0],formulation.ndim,self.number_of_load_increments),dtype=np.float64)
        TotalS = np.zeros((formulation.meshes[2].points.shape[0],formulation.ndim,self.number_of_load_increments),dtype=np.float64)
        # TotalDisp = np.zeros((mesh.points.shape[0],int(formulation.ndim**2),self.number_of_load_increments),dtype=np.float64)

        # PRE-ASSEMBLY
        if caller != "Solve":
            print('Assembling the system and acquiring neccessary information for the analysis...')
        tAssembly=time()

        # APPLY DIRICHELT BOUNDARY CONDITIONS AND GET DIRICHLET RELATED FORCES
        boundary_condition.GetDirichletBoundaryConditions(formulation, mesh, material, solver, self)

        # GET BOUNDARY INFO
        self.GetBoundaryInfo(mesh,formulation,boundary_condition)

        # ALLOCATE FOR GEOMETRY - GetDirichletBoundaryConditions CHANGES THE MESH
        # SO EULERX SHOULD BE ALLOCATED AFTERWARDS
        Eulerx = np.copy(formulation.meshes[0].points)
        Eulerw = np.zeros_like(formulation.meshes[1].points)
        Eulers = np.zeros_like(formulation.meshes[1].points)
        Eulerp = np.zeros((formulation.meshes[0].points.shape[0]))

        # FIND PURE NEUMANN (EXTERNAL) NODAL FORCE VECTOR
        NeumannForces = boundary_condition.ComputeNeumannForces(mesh, material, function_spaces,
            compute_traction_forces=True, compute_body_forces=self.add_self_weight)
        # NeumannForces = np.zeros((mesh.points.shape[0]*formulation.ndim))

        # ASSEMBLE STIFFNESS MATRIX AND TRACTION FORCES
        if self.analysis_type != 'static':
            if formulation.fields == "couple_stress" or formulation.fields == "flexoelectric":
                K, TractionForces, _, M = formulation.Assemble(self, material, Eulerx, Eulerw, Eulers, Eulerp)
            elif formulation.fields == "mechanics" or formulation.fields == "electro_mechanics":
                # STANDARD MECHANINCS ELECTROMECHANICS DYNAMIC ANALYSIS ARE DISPATCHED HERE
                from Florence.FiniteElements.Assembly import Assemble
                fspace = function_spaces[0] if (mesh.element_type=="hex" or mesh.element_type=="quad") else function_spaces[1]
                K, TractionForces, _, M = Assemble(self, fspace, formulation, mesh, material,
                    Eulerx, Eulerp)
        else:
            K, TractionForces = formulation.Assemble(self, material, Eulerx, Eulerw, Eulers, Eulerp)[:2]

        print('Finished all pre-processing stage. Time elapsed was', time()-tAssembly, 'seconds')

        if self.analysis_type != 'static':
            boundary_condition.ConvertStaticsToDynamics(mesh, self.number_of_load_increments)
            TotalDisp, TotalW, TotalS = self.DynamicSolver(formulation, solver,
                K, M, NeumannForces, NodalForces, Residual,
                mesh, TotalDisp, TotalW, TotalS, Eulerx, Eulerw, Eulers, Eulerp, material, boundary_condition)
        else:
            TotalDisp, TotalW, TotalS = self.StaticSolver(formulation, solver,
                K, NeumannForces, NodalForces, Residual,
                mesh, TotalDisp, TotalW, TotalS, Eulerx, Eulerw, Eulers, Eulerp, material, boundary_condition)

        solution = self.__makeoutput__(mesh, TotalDisp, formulation, function_spaces, material)
        if formulation.fields == "couple_stress" or formulation.fields == "flexoelectric":
            solution.solW = TotalW
            solution.solS = TotalS
        return solution

    def StaticSolver(self, formulation, solver, K,
            NeumannForces, NodalForces, Residual,
            mesh, TotalDisp, TotalW, TotalS, Eulerx, Eulerw, Eulers, Eulerp, material, boundary_condition):

        LoadIncrement = self.number_of_load_increments
        LoadFactor = 1./LoadIncrement
        AppliedDirichletInc = np.zeros(boundary_condition.applied_dirichlet.shape[0],dtype=np.float64)
        Increment = 0

        post_process = PostProcess(formulation.ndim,formulation.nvar)
        post_process.SetAnalysis(analysis_type=self.analysis_type, analysis_nature=self.analysis_nature)

        # APPLY NEUMANN BOUNDARY CONDITIONS
        DeltaF = LoadFactor*NeumannForces
        Residual[:,:] += DeltaF

        K_b, F_b, _ = boundary_condition.ApplyDirichletGetReducedMatrices(K,Residual,
                boundary_condition.applied_dirichlet,LoadFactor=LoadFactor)

        for Increment in range(LoadIncrement):
            t_increment=time()
            # SOLVE THE SYSTEM
            sol = solver.Solve(K_b, F_b, reuse_factorisation=True)

            AppliedDirichletInc = LoadFactor*boundary_condition.applied_dirichlet
            dU = post_process.TotalComponentSol(sol, boundary_condition.columns_in,
                boundary_condition.columns_out, AppliedDirichletInc,0,K.shape[0])

            # STORE TOTAL SOLUTION DATA
            TotalDisp[:,:,Increment] += dU

            # REDUCED ACCUMULATED FORCE
            Residual[:,:] += DeltaF
            F_b = boundary_condition.ApplyDirichletGetReducedMatrices(K,Residual,
                AppliedDirichletInc,LoadFactor=1.0,
                only_residual=True)[boundary_condition.columns_in,0]
            # F_b = boundary_condition.GetReducedVectors(Residual,only_residual=True)[0]

            # UPDATE
            Eulerx += dU[:,:formulation.ndim]
            Eulerp += dU[:,-1]
            TotalW[:,:,Increment], TotalS[:,:,Increment] = formulation.GetAugmentedSolution(self,
                material, TotalDisp, Eulerx, Eulerw, Eulers, Eulerp)

            # LOG REQUESTS
            self.LogSave(formulation, TotalDisp, Increment)

            # BREAK AT A SPECIFICED LOAD INCREMENT IF ASKED FOR
            if self.break_at_increment != -1 and self.break_at_increment is not None:
                if self.break_at_increment == Increment:
                    if self.break_at_increment < LoadIncrement - 1:
                        print("\nStopping at increment {} as specified\n\n".format(Increment))
                        TotalDisp = TotalDisp[:,:,:Increment]
                        self.number_of_load_increments = Increment
                    break

            print('Finished Load increment', Increment, 'in', time()-t_increment, 'seconds\n')

        # ADD EACH INCREMENTAL CONTRIBUTION TO MAKE IT CONSISTENT WITH THE NONLINEAR ANALYSIS
        for i in range(TotalDisp.shape[2]-1,0,-1):
            TotalDisp[:,:,i] = np.sum(TotalDisp[:,:,:i+1],axis=2)

        # COMPUTE DISSIPATION OF ENERGY THROUGH TIME
        if self.compute_energy:
            energy_info = self.ComputeEnergy(formulation.function_spaces[0],mesh,material,formulation,Eulerx,Eulerp)
            formulation.strain_energy = energy_info[0]
            formulation.electrical_energy = energy_info[1]

        solver.CleanUp()

        return TotalDisp, TotalW, TotalS

    def DynamicSolver(self, formulation, solver, K, M,
            NeumannForces, NodalForces, Residual,
            mesh, TotalDisp, TotalW, TotalS, Eulerx, Eulerw, Eulers, Eulerp, material, boundary_condition):

        LoadIncrement = self.number_of_load_increments
        LoadFactor = self.total_time/LoadIncrement

        post_process = PostProcess(formulation.ndim,formulation.nvar)
        post_process.SetAnalysis(analysis_type=self.analysis_type, analysis_nature=self.analysis_nature)

        if NeumannForces.ndim == 2 and NeumannForces.shape[1]==1:
            tmp = np.zeros((NeumannForces.shape[0],LoadIncrement))
            tmp[:,0] = NeumannForces[:,0]
            NeumannForces = tmp

        dumU = np.zeros((mesh.points.shape[0]*formulation.nvar))
        dumU[boundary_condition.columns_out] = boundary_condition.applied_dirichlet[:,0]
        TotalDisp[:,:formulation.nvar,0] = dumU.reshape(mesh.points.shape[0],formulation.nvar)
        # INITIALISE VELOCITY AND ACCELERATION
        velocities     = np.zeros((mesh.points.shape[0]*formulation.ndim))
        accelerations  = np.zeros((mesh.points.shape[0]*formulation.ndim))
        # COMPUTE DAMPING MATRIX BASED ON MASS
        D = 0.0
        if self.include_physical_damping:
            D = self.damping_factor*M

        if formulation.fields == "electro_mechanics" or formulation.fields == "flexoelectric":
            # self.GetBoundaryInfo(mesh, formulation,boundary_condition)
            M_mech = M[self.mechanical_dofs,:][:,self.mechanical_dofs]
            if self.include_physical_damping:
                D_mech = D[self.mechanical_dofs,:][:,self.mechanical_dofs]
        else:
            M_mech = M
            D_mech = D

        # COMPUTE INITIAL ACCELERATION FOR TIME STEP 0
        Residual = np.zeros_like(Residual)
        InitResidual = Residual + NeumannForces[:,0][:,None]
        if formulation.fields == "electro_mechanics" or formulation.fields == "flexoelectric":
            accelerations[:] = solver.Solve(M_mech, -InitResidual[self.mechanical_dofs].ravel())
        else:
            accelerations[:] = solver.Solve(M, InitResidual.ravel() )

        # COMPUTE AUGMENTED K (INCLUDES INERTIA EFFECT)
        K          += (self.gamma/self.beta/LoadFactor)*D + (1./self.beta/LoadFactor**2)*M
        # GET REDUCED VARIABLES
        K_b, F_b, _ = boundary_condition.GetReducedMatrices(K,Residual)

        if self.lump_rhs:
            M_mech = M_mech.sum(axis=1).A.ravel() # FOR CSR
            # M_mech = M_mech.sum(axis=0).ravel() # FOR CSC
            if self.include_physical_damping:
                D_mech = D_mech.sum(axis=1).A.ravel()

        for Increment in range(1,LoadIncrement):
            t_increment=time()

            # FIXED INCREMENTAL DIRICHLET
            AppliedDirichletInc = boundary_condition.applied_dirichlet[:,Increment-1]

            # APPLY NEUMANN BOUNDARY CONDITIONS
            DeltaF = NeumannForces[:,Increment][:,None]
            NodalForces = DeltaF

            # ACCUMULATED FORCE
            if self.include_physical_damping:
                if self.lump_rhs:
                    Residual[self.mechanical_dofs,0] = (1./self.beta/LoadFactor**2)*M_mech*TotalDisp[:,:formulation.ndim,Increment-1].ravel() +\
                        (1./self.beta/LoadFactor)*M_mech*velocities + (0.5/self.beta - 1.)*M_mech*accelerations +\
                        (self.gamma/self.beta/LoadFactor)*D_mech*TotalDisp[:,:formulation.ndim,Increment-1].ravel() +\
                        (self.gamma/self.beta - 1.)*D_mech*velocities -\
                        LoadFactor*((1-self.gamma)-self.gamma*(0.5/self.beta - 1.))*D_mech*accelerations
                else:
                    Residual[self.mechanical_dofs,0] = (1./self.beta/LoadFactor**2)*M_mech.dot(TotalDisp[:,:formulation.ndim,Increment-1].ravel()) +\
                        (1./self.beta/LoadFactor)*M_mech.dot(velocities) + (0.5/self.beta - 1.)*M_mech.dot(accelerations) +\
                        (self.gamma/self.beta/LoadFactor)*D_mech.dot(TotalDisp[:,:formulation.ndim,Increment-1].ravel()) +\
                        (self.gamma/self.beta - 1.)*D_mech.dot(velocities) -\
                        LoadFactor*((1-self.gamma)-self.gamma*(0.5/self.beta - 1.))*D_mech.dot(accelerations)
            else:
                if self.lump_rhs:
                    Residual[self.mechanical_dofs,0] = (1./self.beta/LoadFactor**2)*M_mech*TotalDisp[:,:formulation.ndim,Increment-1].ravel() +\
                        (1./self.beta/LoadFactor)*M_mech*velocities + (0.5/self.beta - 1.)*M_mech*accelerations
                else:
                    Residual[self.mechanical_dofs,0] = (1./self.beta/LoadFactor**2)*M_mech.dot(TotalDisp[:,:formulation.ndim,Increment-1].ravel()) +\
                        (1./self.beta/LoadFactor)*M_mech.dot(velocities) + (0.5/self.beta - 1.)*M_mech.dot(accelerations)
            Residual += DeltaF

            # CHECK CONTACT AND ASSEMBLE IF DETECTED
            if self.has_contact:
                Eulerx = mesh.points + TotalDisp[:,:formulation.ndim,Increment-1]
                TractionForcesContact = np.zeros_like(Residual)
                TractionForcesContact = self.contact_formulation.AssembleTractions(mesh,material,Eulerx).ravel()*LoadFactor

                if formulation.fields == "electro_mechanics" or formulation.fields == "flexoelectric":
                    Residual[self.mechanical_dofs,0] -= TractionForcesContact
                elif formulation.fields == "mechanics" or formulation.fields == "couple_stress":
                    Residual[:,0] -= TractionForcesContact
                else:
                    raise NotImplementedError("Contact algorithm for {} is not available".format(formulation.fields))

            # REDUCED ACCUMULATED FORCE
            F_b = boundary_condition.ApplyDirichletGetReducedMatrices(K,Residual,
                boundary_condition.applied_dirichlet[:,Increment],LoadFactor=1.0,
                mass=M,only_residual=True)[boundary_condition.columns_in,0]

            # SOLVE THE SYSTEM
            sol = solver.Solve(K_b, F_b, reuse_factorisation=True)

            dU = post_process.TotalComponentSol(sol, boundary_condition.columns_in,
                boundary_condition.columns_out, AppliedDirichletInc,0,K.shape[0])

            # STORE TOTAL SOLUTION DATA
            TotalDisp[:,:,Increment] += dU

            # UPDATE VELOCITY AND ACCELERATION
            accelerations_old = np.copy(accelerations)
            accelerations = (1./self.beta/LoadFactor**2)*(TotalDisp[:,:formulation.ndim,Increment] -\
                TotalDisp[:,:formulation.ndim,Increment-1]).ravel() -\
                1./self.beta/LoadFactor*velocities + (1.-0.5/self.beta)*accelerations_old
            velocities += LoadFactor*(self.gamma*accelerations + (1-self.gamma)*accelerations_old)

            # UPDATE
            Eulerx += dU[:,:formulation.ndim]
            Eulerp += dU[:,-1]
            if formulation.fields == "couple_stress" or formulation.fields == "flexoelectric":
                TotalW[:,:,Increment], TotalW[:,:,Increment] = formulation.GetAugmentedSolution(self, material,
                    TotalDisp, Eulerx, Eulerw, Eulers, Eulerp)

            # LOG REQUESTS
            self.LogSave(formulation, TotalDisp, Increment)

            # BREAK AT A SPECIFICED LOAD INCREMENT IF ASKED FOR
            if self.break_at_increment != -1 and self.break_at_increment is not None:
                if self.break_at_increment == Increment:
                    if self.break_at_increment < LoadIncrement - 1:
                        print("\nStopping at increment {} as specified\n\n".format(Increment))
                        TotalDisp = TotalDisp[:,:,:Increment]
                        self.number_of_load_increments = Increment
                    break

            # STORE THE INFORMATION IF THE SOLVER BLOWS UP
            if Increment > 0:
                U0 = TotalDisp[:,:,Increment-1].ravel()
                U = TotalDisp[:,:,Increment].ravel()
                tol = 1e200 if Increment < 5 else 10.
                if np.isnan(norm(U)) or np.abs(U.max()/(U0.max()+1e-14)) > tol:
                    print("Solver blew up! Norm of incremental solution is too large")
                    TotalDisp = TotalDisp[:,:,:Increment]
                    self.number_of_load_increments = Increment
                    break

            print('Finished Load increment', Increment, 'in', time()-t_increment, 'seconds\n')

        # COMPUTE DISSIPATION OF ENERGY THROUGH TIME
        if self.compute_energy:
            energy_info = self.ComputeEnergy(formulation.function_spaces[0],mesh,material,formulation,Eulerx,Eulerp)
            formulation.strain_energy = energy_info[0]
            formulation.electrical_energy = energy_info[1]

        solver.CleanUp()

        return TotalDisp, TotalW, TotalS

    def ComputeEnergy(self,function_space,mesh,material,formulation,Eulerx,Eulerp):

        strain_energy = 0.
        electrical_energy = 0.
        for elem in range(mesh.nelem):
            LagrangeElemCoords = mesh.points[mesh.elements[elem,:],:]
            EulerElemCoords = Eulerx[mesh.elements[elem,:],:]
            ElectricPotentialElem = Eulerp[mesh.elements[elem,:]]

            energy = formulation.GetEnergy(function_space, material,
                LagrangeElemCoords, EulerElemCoords, ElectricPotentialElem, self, elem)
            strain_energy += energy[0]
            electrical_energy += energy[1]

        return strain_energy, electrical_energy

    def ComputeEnergyDissipation(self,function_space,mesh,material,formulation,fem_solver,
        Eulerx, TotalDisp, NeumannForces, M, velocities, Increment):

        internal_energy = 0.
        for elem in range(mesh.nelem):
            LagrangeElemCoords = mesh.points[mesh.elements[elem,:],:]
            EulerElemCoords = Eulerx[mesh.elements[elem,:],:]

            internal_energy += formulation.GetEnergy(function_space, material,
                LagrangeElemCoords, EulerElemCoords, fem_solver, elem)

        if formulation.fields == "flexoelectric":
            M_mech = M[self.mechanical_dofs,:][:,self.mechanical_dofs]
            kinetic_energy = 0.5*np.dot(velocities.ravel(),M_mech.dot(velocities.ravel()))
        else:
            kinetic_energy = 0.5*np.dot(velocities.ravel(),M.dot(velocities.ravel()))

        external_energy = np.dot(TotalDisp[:,:,Increment].ravel(),NeumannForces[:,Increment])

        total_energy = internal_energy + kinetic_energy - external_energy
        return total_energy, internal_energy, kinetic_energy, external_energy

    def ComputePowerDissipation(self,function_space,mesh,material,formulation,fem_solver, Eulerx, TotalDisp,
        NeumannForces, M, velocities, accelerations, Increment):

        internal_energy = 0.
        for elem in range(mesh.nelem):
            LagrangeElemCoords = mesh.points[mesh.elements[elem,:],:]
            EulerElemCoords    = Eulerx[mesh.elements[elem,:],:]
            VelocityElem       = velocities[mesh.elements[elem,:],:]

            internal_energy += formulation.GetLinearMomentum(function_space, material,
                LagrangeElemCoords, EulerElemCoords, VelocityElem, fem_solver, elem)

        if formulation.fields == "electro_mechanics":
            M_mech = M[self.mechanical_dofs,:][:,self.mechanical_dofs]
            kinetic_energy = np.dot(velocities.ravel(),M_mech.dot(accelerations.ravel()))
        else:
            kinetic_energy = np.dot(velocities.ravel(),M.dot(accelerations.ravel()))

        external_energy = np.dot(velocities.ravel(),NeumannForces[:,Increment])

        total_energy = internal_energy + kinetic_energy - external_energy
        return total_energy, internal_energy, kinetic_energy, external_energy

max_stars_repo_name: jdlaubrie/florence

max_issues_repo_name: jdlaubrie/florence

max_forks_repo_name: jdlaubrie/florence

reason: import numpy,from numpy,import scipy


text:
from typing import Union

import numpy as np

from fmodeling.seismic.dynamic.zoeppritz_coeffs import pdownpup
from fmodeling.seismic.dynamic.zoeppritz_coeffs_water import pdownpup_water
from objects.seismic.boundaries.base_boundary import BaseBoundary
from objects.seismic.boundaries.base_boundary import register_boundary
from objects.seismic.waves import WD, WT

class ReflectionUp(BaseBoundary):
    @property
    def wave_direction_1(self):
        return WD.DOWN

    @property
    def wave_direction_2(self):
        return WD.UP

class PDownPUpReflection(ReflectionUp):
    @property
    def wave_type_1(self):
        return WT.P

    @property
    def wave_type_2(self):
        return WT.P

    @classmethod
    def calculate_coeff(cls, vp1: Union[float, np.ndarray], vs1: Union[float, np.ndarray], rho1: Union[float, np.ndarray],
                            vp2: Union[float, np.ndarray], vs2: Union[float, np.ndarray], rho2: Union[float, np.ndarray],
                        theta1: Union[float, np.ndarray] = 0):
        return pdownpup(vp1, vs1, rho1, vp2, vs2, rho2, theta1)

@register_boundary
class PDownPUpWaterReflection(ReflectionUp):
    @property
    def wave_type_1(self):
        return WT.P

    @property
    def wave_type_2(self):
        return WT.P

    @classmethod
    def calculate_coeff(cls, vp1: Union[float, np.ndarray], vs1: Union[float, np.ndarray],
                        rho1: Union[float, np.ndarray],
                        vp2: Union[float, np.ndarray], vs2: Union[float, np.ndarray], rho2: Union[float, np.ndarray],
                        theta1: Union[float, np.ndarray] = 0):
        return pdownpup_water(vp1, rho1, vp2, vs2, rho2, theta1)

ext: py

max_stars_repo_name: lomonosouv/rpsi

max_issues_repo_name: lomonosouv/rpsi

max_forks_repo_name: lomonosouv/rpsi

reason: import numpy


text:
#Copyright (C) 2020. Huawei Technologies Co., Ltd. All rights reserved.

#This program is free software; you can redistribute it and/or modify it under the terms of the BSD 3-Clause License.

#This program is distributed in the hope that it will be useful, but WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the BSD 3-Clause License for more details.

import os
import sys
import argparse
from datetime import datetime

import numpy as np
import torch
import torch.nn as nn
import torch.optim as optim
import torchvision
import torchvision.transforms as transforms

from torch.utils.data import DataLoader

from torch.autograd import Variable
import torch.nn.functional as F
from models.resnet import ResNet_disout

from utils import get_training_dataloader, get_test_dataloader, WarmUpLR

parser = argparse.ArgumentParser()

parser.add_argument('--data_root', type=str,default='/cache/cifar10/')
parser.add_argument('-w', type=int, default=2, help='number of workers for dataloader')
parser.add_argument('-b', type=int, default=128, help='batch size for dataloader')
parser.add_argument('-s', type=bool, default=True, help='whether shuffle the dataset')
parser.add_argument('-warm', type=int, default=1, help='warm up training phase')
parser.add_argument('-lr', type=float, default=0.1, help='initial learning rate')
parser.add_argument('--epochs', type=int, default=200)
parser.add_argument('--dist_prob', type=float, default=0.09)
parser.add_argument('--block_size', type=int, default=6)
parser.add_argument('--alpha', type=float, default=5)
args,unparsed = parser.parse_known_args()

def train(epoch):

    train_loss = 0.0
    correct = 0.0

    net.train()
    for batch_index, (images, labels) in enumerate(training_loader):
        if epoch <= args.warm:
            warmup_scheduler.step()

        images = Variable(images)
        labels = Variable(labels)

        labels = labels.cuda()
        images = images.cuda()

        batch_size=images.size(0)

        optimizer.zero_grad()
        outputs = net(images)

        loss = loss_function(outputs, labels)
        loss.backward()
        optimizer.step()

        train_loss += loss.item()
        _, preds = outputs.max(1)
        correct += preds.eq(labels).sum()

        n_iter = (epoch - 1) * len(training_loader) + batch_index + 1

    print('Train set Epoch: {epoch}  Average_loss: {loss:.4f}, Accuracy: {acc:.4f}'.format(epoch=epoch,
        loss=train_loss / len(training_loader.dataset),
        acc=correct.float() / len(training_loader.dataset)))

def eval(epoch):
    net.eval()

    test_loss = 0.0
    correct = 0.0

    for (images, labels) in cifar10_test_loader:
        images = Variable(images)
        labels = Variable(labels)

        images = images.cuda()
        labels = labels.cuda()

        outputs = net(images)
        loss = loss_function(outputs, labels)
        test_loss += loss.item()
        _, preds = outputs.max(1)
        correct += preds.eq(labels).sum()

    print('Test set: Average_loss: {:.4f}, Accuracy: {:.4f}'.format(
        test_loss / len(cifar10_test_loader.dataset),
        correct.float() / len(cifar10_test_loader.dataset)
    ))

    return correct.float() / len(cifar10_test_loader.dataset)

if __name__ == '__main__':

    training_loader = get_training_dataloader(args.data_root,
        (0.4914, 0.4822, 0.4465),
        (0.2023, 0.1994, 0.2010),
        num_workers=args.w,
        batch_size=args.b,
        shuffle=args.s
    )

    cifar10_test_loader = get_test_dataloader(args.data_root,
        (0.4914, 0.4822, 0.4465),
        (0.2023, 0.1994, 0.2010),
        num_workers=args.w,
        batch_size=args.b,
        shuffle=args.s
    )

    net = ResNet_disout(depth=56, num_classes=10,dist_prob=args.dist_prob,block_size=args.block_size,alpha=args.alpha,
                        nr_steps=len(training_loader)*args.epochs).cuda()

    loss_function = nn.CrossEntropyLoss()
    optimizer = optim.SGD(net.parameters(), lr=args.lr, momentum=0.9, weight_decay=5e-4)
    train_scheduler = optim.lr_scheduler.MultiStepLR(optimizer, milestones=[60, 120, 160], gamma=0.2) #learning rate decay
    iter_per_epoch = len(training_loader)
    warmup_scheduler = WarmUpLR(optimizer, iter_per_epoch * args.warm)

    best_acc = 0.0
    for epoch in range(1, args.epochs+1):
        if epoch > args.warm:
            train_scheduler.step(epoch)
        train(epoch)
        acc = eval(epoch)
        if best_acc < acc:
            best_acc=acc
    print('best_acc:',best_acc)

max_stars_repo_name: huawei-noah/Disout

max_issues_repo_name: huawei-noah/Disout

max_forks_repo_name: huawei-noah/Disout


text:
################################################################################
######################### CREATING VARIABLE CLASS BELOW ########################
###### No need to edit below unless making major changes to code ###############
from config import *
import os
import sys
from datetime import datetime
import numpy as np

# Create objection containing all the variables set above
class create_config_vars:
    def __init__(self,**kwargs):

        self.output_dir = kwargs.get('output_dir','')
        self.satellite_product = kwargs.get('satellite_product','')
        self.satellite_file_dir = kwargs.get('satellite_file_dir','')
        self.satellite_file_suffix = kwargs.get('satellite_file_suffix','nc')
        self.apply_quality_control = kwargs.get('apply_quality_control',False)
        self.lat_var_name = kwargs.get('lat_var_name','')
        self.lon_var_name = kwargs.get('lon_var_name','')
        self.lev_var_name = kwargs.get('lev_var_name','')
        self.ak_lev_var_name = kwargs.get('ak_lev_var_name',self.lev_var_name)
        self.time_var_name = kwargs.get('time_var_name','')
        self.ak_var_name = kwargs.get('ak_var_name','')
        self.prior_var_name = kwargs.get('prior_var_name','')
        self.o3_var_name = kwargs.get('o3_var_name','')
        self.model_repo_url = kwargs.get('model_repo_url','')
        self.keep_model_downloads = kwargs.get('keep_model_downloads',True)
        self.verbose = kwargs.get('verbose',True)
        self.model_temporal_res = kwargs.get('model_temporal_res', '')
        self.start_date = datetime(*kwargs.get('start_date',(2000,1,1)))
        self.end_date = datetime(*kwargs.get('end_date',(2020,1,1)))

        model_download_path = os.path.join(self.output_dir, 'model_download')
        self.model_download_path = model_download_path

# Make sure the model url ends with a slash otherwise the URL join doesn't work
if model_repo_url[-1] != '/':
    model_repo_url = model_repo_url + '/'
# Make sure the model temporal resolution makes sense
try:
    num_res = int(model_temporal_res.split()[0])
except (NameError, SyntaxError) as e:
    print("** Model temporal resoltion ({}) not a valid input./n Format should be N hourly, N daily or N monthly".format(
        model_temporal_res
    ))
    sys.exit()
if model_temporal_res.split()[1].lower() not in ['hourly','daily','monthly']:
    print("** Model temporal resoltion ({}) not a valid input./n Format should be N hourly, N daily or N monthly".format(
        model_temporal_res
    ))
    sys.exit()

if np.sum(list(satellite_product.values())) > 1:
    print("** More than one satellite product has been selected in config.py. Only choose one.")
    sys.exit()

elif np.sum(list(satellite_product.values())) < 1:
    print("** No satellite product has been selected. Choose one.")
    sys.exit()

else:
    satellite_product_name = [k for k,v in satellite_product.items() if v][0]

print("### SAMPLING MODEL FOR {} ###".format(satellite_product_name))

if ak_level_var_name == None:
    ak_level_var_name = level_var_name

config_vars = create_config_vars(output_dir = output_dir,
                                satellite_product = satellite_product_name,
                                satellite_file_dir =  satellite_file_dir,
                                start_date = start_date,
                                end_date = end_date,
                                satellite_file_suffix = satellite_file_suffix,
                                apply_quality_control = apply_quality_control,
                                lat_var_name = latitude_var_name,
                                lon_var_name = longitude_var_name,
                                lev_var_name = level_var_name,
                                ak_lev_var_name = ak_level_var_name,
                                time_var_name = time_var_name,
                                o3_var_name = o3_var_name,
                                ak_var_name = ak_var_name,
                                prior_var_name = prior_var_name,
                                model_repo_url = model_repo_url,
                                keep_model_downloads = keep_model_downloads,
                                verbose = verbose,
                                model_temporal_res = model_temporal_res)


text:
# from __future__ import absolute_import
# from __future__ import division
# from __future__ import print_function

import glob
import os.path as osp
import numpy as np
# import json_tricks as json
import json
import pickle
import logging
import os
import copy
from torch.utils.data import Dataset
import cv2
import torch
import torchvision.transforms as transforms
import glob
from tqdm import tqdm
import numpy as np
import sys
if sys.version_info[0] == 2:
    import xml.etree.cElementTree as ET
else:
    import xml.etree.ElementTree as ET
namestr = ['CVLab1','CVLab2','CVLab3','CVLab4','IDIAP1','IDIAP2','IDIAP3']

CONNS =  [[0, 1],
         [0, 2],
         [0, 3],
         [3, 4],
         [4, 5],
         [0, 9],
         [9, 10],
         [10, 11],
         [2, 6],
         [2, 12],
         [6, 7],
         [7, 8],
         [12, 13],
         [13, 14]] # only 15 keypoints

class Wildtrack:
    def __init__(self, image_folder, keypoint_folder, view_set,is_train = True): # TODO add keypoint foldercfg,
        self.view_set = view_set
        # self.cam_list = [(0,x) for x in self.view_set] # get the HD images (0,)
        self.view_num = len(self.view_set)
        self.image_folder = image_folder
        # self.transform = transform
        self.conns = CONNS
        self.image_size = np.array([1920,1080])
        # self.input_size = cfg.dataset.input_size
        # self.input_size = np.array([960,512])
        # self.heatmap_size = self.input_size / 2
        # self.heatmap_size = self.heatmap_size.astype(np.int16)

        self.num_joints = 15 #cfg.NETWORK.NUM_JOINTS
        self.paf_num = len(CONNS)

        # self.sigma = 4 #cfg.NETWORK.SIGMA
        # self.single_size = 512*424
        self.istrain = is_train
        # 读取k_calibration, ksync, 以depth 为准对齐一次即可
        # self.scene_num = len(self.scene_list)

        # calculate the total frame idx for the specific idx
        # read in the keypoint file as the order
        self.kpfiles = os.listdir(keypoint_folder)
        self.kpfiles.sort(key=lambda x:int(x.split('.')[0]))

        self.camera_parameters = self._get_cam()
        self.max_people = 40

    def _get_cam(self):
        cameras = {}
        for cam in self.view_set:
            sel_cam = {}
            root = ET.parse('/home/fbh/wildtrackdataset/Wildtrack_dataset/calibrations/extrinsic/extr_'+namestr[cam]+'.xml').getroot()
            R = root[0].text
            R = R.replace('\n', '')
            R = R.split(' ')
            Rm = []
            for str in R:
                if str!='':
                    Rm.append(float(str))
            rx = Rm[0]
            Rx = np.zeros([3,3])
            Rx[0,0] = 1
            Rx[1,1] = np.cos(rx)
            Rx[2,1] = np.sin(rx)
            Rx[1,2] = -np.sin(rx)
            Rx[2,2] = np.cos(rx)

            ry = Rm[1]
            Ry = np.zeros([3, 3])
            Ry[0, 0] = np.cos(ry)
            Ry[1, 1] = 1
            Ry[2, 0] = -np.sin(ry)
            Ry[0, 2] = np.sin(ry)
            Ry[2, 2] = np.cos(ry)

            rz = Rm[2]
            Rz = np.zeros([3, 3])
            Rz[0, 0] = np.cos(rz)
            Rz[0, 1] = -np.sin(rz)
            Rz[1, 0] = np.sin(rz)
            Rz[1, 1] = np.cos(rz)
            Rz[2, 2] = 1
            Rs = np.dot(np.dot(Rx,Ry),Rz)
            sel_cam['R'] = Rs
            # extr[:,0:3,i] = Rs
            t = root[1].text
            t = t.replace('\n', '')
            t = t.split(' ')
            tm = []
            for str in t:
                if str!='':
                    tm.append(float(str))
            Ts = np.stack(tm)
            sel_cam['t'] = Ts.reshape((3, 1))
            # extr[:,3,i] = np.stack(tm)
            sel_cam['distCoef'] = np.zeros(5)
            root = ET.parse('/home/fbh/wildtrackdataset/Wildtrack_dataset/calibrations/intrinsic_zero/intr_' + namestr[cam] + '.xml').getroot() # zero
            R = root[0][3].text
            R = R.replace('\n', '')
            R = R.split(' ')
            Rm = []
            for str in R:
                if str != '':
                    Rm.append(float(str))
            Rm = np.stack(Rm,axis = 0)
            Rm = np.reshape(Rm,[3,3])
            sel_cam['K'] = Rm
            # intr[:,:,i] = Rm

    def __len__(self):
        return int(np.sum(self.num_pers))

    def __generate_meta__(self):  # get the transfer (local to global) and the distosion parameters
        meta = dict()
        meta['root'] = list()

        for index in tqdm(range(self.dataset_len)):
            findtable = self.until_sum - (index + 1)
            scene_index = np.min(np.where(findtable>=0)[0])
            until_num = np.sum(self.num_pers[:scene_index])
            kp_idx = index - until_num
            kp_file = self.anno_kp_files[scene_index][kp_idx]

            with open(kp_file,'rb') as kp:
                try:
                    kp_row_data = json.load(kp)
                except:
                    continue
            kp3d_body_data = kp_row_data['bodies']
            nposes = len(kp3d_body_data)
            if nposes == 0:
                continue
            number_people = torch.tensor(nposes)
            camera_paras = self.camera_parameters[scene_index]

            per_info_media = dict()
            for cam_node in self.cam_list:
                # generate meta file
                # different views

                per_info = dict()
                per_info['img_height'] = int(self.image_size[1])
                per_info['img_width'] = int(self.image_size[0])
                per_info['dataset'] = 'CMUP'
                if self.istrain:
                    per_info['isValidation'] = 0
                else:
                    per_info['isValidation'] = 1

                view_para = camera_paras[cam_node]
                K = view_para['K']
                R = view_para['R']
                T = view_para['t']
                Kd = view_para['distCoef'] # output this Kd for recovering
                prefix = '{:02d}_{:02d}'.format(cam_node[0], cam_node[1])
                postfix = osp.basename(kp_file).replace('body3DScene', '')
                file_name = osp.join(self.image_folder,self.scene_list[scene_index],'hdImgs',prefix,
                                        prefix+postfix)
                file_name = file_name.replace('json','jpg')
                per_info['img_paths'] = file_name

                pose_3d = []
                pose_info = []
                for n in range(nposes):
                    pose3d = np.array(kp3d_body_data[n]['joints19']).reshape((-1, 4))
                    # process the joint into 15 keypoints
                    pose3d = pose3d[:15,:].copy() # only consider 15 keypoints # only 15 keypoints
                    anno_vis = 2 * (pose3d[:, -1] > 0.1)
                    pose3d_proc = pose3d[...,:3].copy()
                    pose2d, depth_val, points_3d_cam, fx, fy, cx, cy = self.__projectjointsPoints__(pose3d_proc.transpose(),K, R, T, Kd) # get the corresponding joint depth value

                    x_check = np.bitwise_and(pose2d[:, 0] >= 0,
                                                pose2d[:, 0] <= self.image_size[0] - 1) #(15,) bool
                    y_check = np.bitwise_and(pose2d[:, 1] >= 0,
                                                pose2d[:, 1] <= self.image_size[1] - 1)

                    check = np.bitwise_and(x_check, y_check) # check bool se
                    anno_vis[np.logical_not(check)] = 0

                    # # process the joints
                    for i in range(len(pose2d)):
                        if (np.min(pose2d[i, :2]) < 0 or
                                pose2d[i, 0] >= self.image_size[0] or
                                pose2d[i, 1] >= self.image_size[1]):
                            anno_vis[i] = 0
                    anno_vis = np.expand_dims(anno_vis,axis = -1)
                    joints_info = np.concatenate([pose2d, depth_val, anno_vis, points_3d_cam, fx, fy, cx, cy], axis=-1) # the joint info
                    pose_info.append(joints_info[None,...])
                    pose_3d.append(pose3d[None,...])

                pose_3d = np.concatenate(pose_3d, axis=0)
                per_info['bodys_3d'] = pose_3d
                per_info['cam'] = view_para # contains K R T Kd
                pose_info = np.concatenate(pose_info,axis=0)
                per_info['bodys'] = pose_info
                per_info_media[cam_node] = per_info
                # meta['root'].append(per_info)
            meta['root'].append(per_info_media)

        if self.istrain:
            writen_file = osp.join(self.image_folder,'cmu_data_train_new5_multi.pkl')
        else:
            writen_file = osp.join(self.image_folder,'cmu_data_test_alterview.pkl')

        with open(writen_file,'wb') as f:
            pickle.dump(meta, f)

    def __projectjointsPoints__(self, X, K, R, t, Kd):
        """
        Projects points X (3xN) using camera intrinsics K (3x3),
        extrinsics (R,t) and distortion parameters Kd=[k1,k2,p1,p2,k3].
        Roughly, x = K*(R*X + t) + distortion
        See http://docs.opencv.org/2.4/doc/tutorials/calib3d/camera_calibration/camera_calibration.html
        or cv2.projectPoints
        """

        x = (np.dot(R, X) + t)  # panoptic to kinect color scaling  cm to m
        points_3d_cam = x.copy().transpose() #[J,3]
        depth_val = x[2:3,:].transpose() # the depth value of current joint, metrix is meter [N,1]

        x[0:2, :] = x[0:2, :] / (x[2, :] + 1e-5)

        r = x[0, :] * x[0, :] + x[1, :] * x[1, :]

        # 去畸变
        x[0, :] = x[0, :] * (1 + Kd[0] * r + Kd[1] * r * r + Kd[4] * r * r * r
                            ) + 2 * Kd[2] * x[0, :] * x[1, :] + Kd[3] * (
                                r + 2 * x[0, :] * x[0, :])
        x[1, :] = x[1, :] * (1 + Kd[0] * r + Kd[1] * r * r + Kd[4] * r * r * r
                            ) + 2 * Kd[3] * x[0, :] * x[1, :] + Kd[2] * (
                                r + 2 * x[1, :] * x[1, :])

        x[0, :] = K[0, 0] * x[0, :] + K[0, 1] * x[1, :] + K[0, 2]
        x[1, :] = K[1, 0] * x[0, :] + K[1, 1] * x[1, :] + K[1, 2]

        # depth_val_norm = depth_val * W / f  # absolute depth sensing
        pose_2d = x[:2,:].copy().transpose()
        fx = K[0,0] * np.ones((pose_2d.shape[0],1))
        fy = K[1,1] * np.ones((pose_2d.shape[0],1))
        cx = K[0,2] * np.ones((pose_2d.shape[0],1))
        cy = K[1,2] * np.ones((pose_2d.shape[0],1))

        return pose_2d, depth_val, points_3d_cam, fx, fy, cx, cy

if __name__ == '__main__':
    img_path = '/Extra/fanbohao/fbh_data/wildtrackdataset/Wildtrack_dataset/Image_subsets/'
    kp_path = '/Extra/fanbohao/fbh_data/wildtrackdataset/Wildtrack_dataset/annotations_positions/'
    view_set = [0,1,2,3,4,5,6]

    depth_data_train = Wildtrack(img_path, kp_path, view_set,is_train = False)
    depth_data_train.__generate_meta__()

max_stars_repo_name: Yu-Yy/MathPoseGNN

max_issues_repo_name: Yu-Yy/MathPoseGNN

max_forks_repo_name: Yu-Yy/MathPoseGNN


text:
"""Take an image and average all pixel values to find brightness. Used for sync and swim.

Attributes:
    cam_l (TYPE): Description
    cam_r (TYPE): Description
    dur (TYPE): Description
"""
import RPi.GPIO as GPIO
GPIO.setwarnings(False)
GPIO.setmode(GPIO.BCM)
import numpy as np
import time
import math

from lib_utils import *
from lib_leds import LEDS
from lib_camera import Camera

def observe():
    flash = False
    thresh_flash = 100

    cam_r.capture()
    avg_bright_r = np.mean(cam_r.img)
    print('avg_bright_r = {}'.format(avg_bright_r))

    cam_l.capture()
    avg_bright_l = np.mean(cam_l.img)
    print('avg_bright_l = {}'.format(avg_bright_l))

    if avg_bright_r > thresh_flash or avg_bright_l > thresh_flash:
        flash = True

    return flash

cam_r = Camera('right')
cam_l = Camera('left')

dur = np.zeros((100,1))

for i in range(100):
	t_start = time.time()
	observe()
	dur[i] = time.time() - t_start

print(np.mean(dur))
print(np.max(dur))

max_stars_repo_name: fberlinger/BlueSwarm

max_issues_repo_name: fberlinger/BlueSwarm

max_forks_repo_name: fberlinger/BlueSwarm


text:
import io
import os
from fractions import Fraction

import pytest

import asdf
from asdf import types
from asdf import extension
from asdf import util
from asdf import versioning
from asdf.exceptions import AsdfWarning, AsdfConversionWarning

from . import helpers, CustomTestType, CustomExtension

TEST_DATA_PATH = str(helpers.get_test_data_path(''))

class Fractional2dCoord:
    def __init__(self, x, y):
        self.x = x
        self.y = y

class FractionWithInverse(Fraction):
    def __init__(self, *args, **kwargs):
        self._inverse = None

    @property
    def inverse(self):
        return self._inverse

    @inverse.setter
    def inverse(self, value):
        self._inverse = value

class FractionWithInverseType(asdf.CustomType):
    name = 'fraction_with_inverse'
    organization = 'nowhere.org'
    version = (1, 0, 0)
    standard = 'custom'
    types = [FractionWithInverse]

    @classmethod
    def to_tree(cls, node, ctx):
        return {
            "numerator": node.numerator,
            "denominator": node.denominator,
            "inverse": node.inverse
        }

    @classmethod
    def from_tree(cls, tree, ctx):
        result = FractionWithInverse(
            tree["numerator"],
            tree["denominator"]
        )
        yield result
        result.inverse = tree["inverse"]

class FractionWithInverseExtension(CustomExtension):
    @property
    def types(self):
        return [FractionWithInverseType]

    @property
    def tag_mapping(self):
        return [('tag:nowhere.org:custom',
                    'http://nowhere.org/schemas/custom{tag_suffix}')]

    @property
    def url_mapping(self):
        return [('http://nowhere.org/schemas/custom/',
                    util.filepath_to_url(TEST_DATA_PATH) + '/{url_suffix}.yaml')]

def fractiontype_factory():

    class FractionType(types.CustomType):
        name = 'fraction'
        organization = 'nowhere.org'
        version = (1, 0, 0)
        standard = 'custom'
        types = [Fraction]
        handle_dynamic_subclasses = True

        @classmethod
        def to_tree(cls, node, ctx):
            return [node.numerator, node.denominator]

        @classmethod
        def from_tree(cls, tree, ctx):
            return Fraction(tree[0], tree[1])

    return FractionType

def fractional2dcoordtype_factory():

    FractionType = fractiontype_factory()

    class Fractional2dCoordType(types.CustomType):
        name = 'fractional_2d_coord'
        organization = 'nowhere.org'
        standard = 'custom'
        version = (1, 0, 0)
        types = [Fractional2dCoord]

        @classmethod
        def to_tree(cls, node, ctx):
            return {
                "x": node.x,
                "y": node.y
            }

        @classmethod
        def from_tree(cls, tree, ctx):
            return Fractional2dCoord(tree["x"], tree["y"])

    class Fractional2dCoordExtension(CustomExtension):
        @property
        def types(self):
            return [FractionType, Fractional2dCoordType]

    return FractionType, Fractional2dCoordType, Fractional2dCoordExtension

def test_custom_tag():

    FractionType = fractiontype_factory()

    class FractionExtension(CustomExtension):
        @property
        def types(self):
            return [FractionType]

    class FractionCallable(FractionExtension):
        @property
        def tag_mapping(self):
            def check(tag):
                prefix = 'tag:nowhere.org:custom'
                if tag.startswith(prefix):
                    return 'http://nowhere.org/schemas/custom' + tag[len(prefix):]
            return [check]

    yaml = """
a: !<tag:nowhere.org:custom/fraction-1.0.0>
  [2, 3]
b: !core/complex-1.0.0
  0j
    """

    buff = helpers.yaml_to_asdf(yaml)
    with asdf.open(buff, extensions=FractionExtension()) as ff:
        assert ff.tree['a'] == Fraction(2, 3)

        buff = io.BytesIO()
        ff.write_to(buff)

    buff = helpers.yaml_to_asdf(yaml)
    with asdf.open(buff, extensions=FractionCallable()) as ff:
        assert ff.tree['a'] == Fraction(2, 3)

        buff = io.BytesIO()
        ff.write_to(buff)
        buff.close()

def test_version_mismatch():
    yaml = """
a: !core/complex-42.0.0
  0j
    """

    buff = helpers.yaml_to_asdf(yaml)
    with pytest.warns(AsdfConversionWarning, match="tag:stsci.edu:asdf/core/complex"):
        with asdf.open(buff, ignore_version_mismatch=False) as ff:
            assert isinstance(ff.tree['a'], complex)

    # Make sure warning is repeatable
    buff.seek(0)
    with pytest.warns(AsdfConversionWarning, match="tag:stsci.edu:asdf/core/complex"):
        with asdf.open(buff, ignore_version_mismatch=False) as ff:
            assert isinstance(ff.tree['a'], complex)

    # Make sure the warning does not occur if it is being ignored (default)
    buff.seek(0)
    with helpers.assert_no_warnings(AsdfConversionWarning):
        with asdf.open(buff) as ff:
            assert isinstance(ff.tree['a'], complex)

    # If the major and minor match, but the patch doesn't, there
    # should still be a warning.
    yaml = """
a: !core/complex-1.0.1
  0j
    """

    buff = helpers.yaml_to_asdf(yaml)
    with pytest.warns(AsdfConversionWarning, match="tag:stsci.edu:asdf/core/complex"):
        with asdf.open(buff, ignore_version_mismatch=False) as ff:
            assert isinstance(ff.tree['a'], complex)

def test_version_mismatch_file(tmpdir):
    testfile = os.path.join(str(tmpdir), 'mismatch.asdf')
    yaml = """
a: !core/complex-42.0.0
  0j
    """

    buff = helpers.yaml_to_asdf(yaml)
    with open(testfile, 'wb') as handle:
        handle.write(buff.read())

    expected_uri = util.filepath_to_url(str(testfile))

    with pytest.warns(AsdfConversionWarning, match="tag:stsci.edu:asdf/core/complex"):
        with asdf.open(testfile, ignore_version_mismatch=False) as ff:
            assert ff._fname == expected_uri
            assert isinstance(ff.tree['a'], complex)

def test_version_mismatch_with_supported_versions():
    """Make sure that defining the supported_versions field eliminates
    the schema mismatch warning."""

    class CustomFlow:
        pass

    class CustomFlowType(CustomTestType):
        version = '1.1.0'
        supported_versions = ['1.0.0', '1.1.0']
        name = 'custom_flow'
        organization = 'nowhere.org'
        standard = 'custom'
        types = [CustomFlow]

    class CustomFlowExtension(CustomExtension):
        @property
        def types(self):
            return [CustomFlowType]

    yaml = """
flow_thing:
  !<tag:nowhere.org:custom/custom_flow-1.0.0>
    c: 100
    d: 3.14
"""
    buff = helpers.yaml_to_asdf(yaml)
    with helpers.assert_no_warnings():
        asdf.open(buff, ignore_version_mismatch=False,
            extensions=CustomFlowExtension())

def test_versioned_writing(monkeypatch):
    from ..tags.core.complex import ComplexType

    # Create a bogus version map
    monkeypatch.setitem(versioning._version_map, '42.0.0', {
        'FILE_FORMAT': '42.0.0',
        'YAML_VERSION': '1.1',
        'tags': {
            'tag:stsci.edu:asdf/core/complex': '42.0.0',
            'tag:stscu.edu:asdf/core/asdf': '1.0.0'
        },
        # We need to insert these explicitly since we're monkeypatching
        'core': {
            'tag:stsci.edu:asdf/core/complex': '42.0.0',
            'tag:stscu.edu:asdf/core/asdf': '1.0.0'
        },
        'standard': {}
    })

    # Add bogus version to supported versions
    monkeypatch.setattr(versioning, 'supported_versions',
        versioning.supported_versions + [versioning.AsdfVersion('42.0.0')]
    )

    class FancyComplexType(types.CustomType):
        name = 'core/complex'
        organization = 'stsci.edu'
        standard = 'asdf'
        version = (42, 0, 0)
        types = [complex]

        @classmethod
        def to_tree(cls, node, ctx):
            return ComplexType.to_tree(node, ctx)

        @classmethod
        def from_tree(cls, tree, ctx):
            return ComplexType.from_tree(tree, ctx)

    class FancyComplexExtension:
        @property
        def types(self):
            return [FancyComplexType]

        @property
        def tag_mapping(self):
            return []

        @property
        def url_mapping(self):
            return [('http://stsci.edu/schemas/asdf/core/complex-42.0.0',
                     util.filepath_to_url(TEST_DATA_PATH) +
                     '/complex-42.0.0.yaml')]

    tree = {'a': complex(0, -1)}

    buff = io.BytesIO()
    ff = asdf.AsdfFile(tree, version="42.0.0",
                       extensions=[FancyComplexExtension()])
    ff.write_to(buff)

    assert b'complex-42.0.0' in buff.getvalue()

def test_longest_match():
    class FancyComplexExtension:
        @property
        def types(self):
            return []

        @property
        def tag_mapping(self):
            return []

        @property
        def url_mapping(self):
            return [('http://stsci.edu/schemas/asdf/core/',
                     'FOOBAR/{url_suffix}')]

    l = extension.AsdfExtensionList(
        [extension.BuiltinExtension(), FancyComplexExtension()])

    assert l.url_mapping(
        'http://stsci.edu/schemas/asdf/core/asdf-1.0.0') == 'FOOBAR/asdf-1.0.0'
    assert l.url_mapping(
        'http://stsci.edu/schemas/asdf/transform/transform-1.0.0') != 'FOOBAR/transform-1.0.0'

def test_module_versioning():
    class NoModuleType(types.CustomType):
        # It seems highly unlikely that this would be a real module
        requires = ['qkjvqdja']

    class HasCorrectPytest(types.CustomType):
        # This means it requires 1.0.0 or greater, so it should succeed
        requires = ['pytest-1.0.0']

    class DoesntHaveCorrectPytest(types.CustomType):
        requires = ['pytest-91984.1.7']

    nmt = NoModuleType()
    hcp = HasCorrectPytest()
    # perhaps an unfortunate acroynm
    dhcp = DoesntHaveCorrectPytest()

    assert nmt.has_required_modules == False
    assert hcp.has_required_modules == True
    assert dhcp.has_required_modules == False

def test_undefined_tag():
    # This tests makes sure that ASDF still returns meaningful structured data
    # even when it encounters a schema tag that it does not specifically
    # implement as an extension
    from numpy import array

    yaml = """
undefined_data:
  !<tag:nowhere.org:custom/undefined_tag-1.0.0>
    - 5
    - {'message': 'there is no tag'}
    - !core/ndarray-1.0.0
      [[1, 2, 3], [4, 5, 6]]
    - !<tag:nowhere.org:custom/also_undefined-1.3.0>
        - !core/ndarray-1.0.0 [[7],[8],[9],[10]]
        - !core/complex-1.0.0 3.14j
"""
    buff = helpers.yaml_to_asdf(yaml)
    with pytest.warns(None) as warning:
        afile = asdf.open(buff)
        missing = afile.tree['undefined_data']

    assert missing[0] == 5
    assert missing[1] == {'message': 'there is no tag'}
    assert (missing[2] == array([[1, 2, 3], [4, 5, 6]])).all()
    assert (missing[3][0] == array([[7],[8],[9],[10]])).all()
    assert missing[3][1] == 3.14j

    # There are two undefined tags, so we expect two warnings
    assert len(warning) == 2
    for i, tag in enumerate(["also_undefined-1.3.0", "undefined_tag-1.0.0"]):
        assert str(warning[i].message) == (
            "tag:nowhere.org:custom/{} is not recognized, converting to raw "
            "Python data structure".format(tag))

    # Make sure no warning occurs if explicitly ignored
    buff.seek(0)
    with helpers.assert_no_warnings():
        afile = asdf.open(buff, ignore_unrecognized_tag=True)

def test_newer_tag():
    # This test simulates a scenario where newer versions of CustomFlow
    # provides different keyword parameters that the older schema and tag class
    # do not account for. We want to test whether ASDF can handle this problem
    # gracefully and still provide meaningful data as output. The test case is
    # fairly contrived but we want to test whether ASDF can handle backwards
    # compatibility even when an explicit tag class for different versions of a
    # schema is not available.
    class CustomFlow:
        def __init__(self, c=None, d=None):
            self.c = c
            self.d = d

    class CustomFlowType(types.CustomType):
        version = '1.1.0'
        name = 'custom_flow'
        organization = 'nowhere.org'
        standard = 'custom'
        types = [CustomFlow]

        @classmethod
        def from_tree(cls, tree, ctx):
            kwargs = {}
            for name in tree:
                kwargs[name] = tree[name]
            return CustomFlow(**kwargs)

        @classmethod
        def to_tree(cls, data, ctx):
            return dict(c=data.c, d=data.d)

    class CustomFlowExtension(CustomExtension):
        @property
        def types(self):
            return [CustomFlowType]

    new_yaml = """
flow_thing:
  !<tag:nowhere.org:custom/custom_flow-1.1.0>
    c: 100
    d: 3.14
"""
    new_buff = helpers.yaml_to_asdf(new_yaml)
    new_data = asdf.open(new_buff, extensions=CustomFlowExtension())
    assert type(new_data.tree['flow_thing']) == CustomFlow

    old_yaml = """
flow_thing:
  !<tag:nowhere.org:custom/custom_flow-1.0.0>
    a: 100
    b: 3.14
"""
    old_buff = helpers.yaml_to_asdf(old_yaml)
    # We expect this warning since it will not be possible to convert version
    # 1.0.0 of CustomFlow to a CustomType (by design, for testing purposes).
    with pytest.warns(AsdfConversionWarning, match="Failed to convert tag:nowhere.org:custom/custom_flow-1.0.0"):
        asdf.open(old_buff, extensions=CustomFlowExtension())

def test_incompatible_version_check():
    class TestType0(types.CustomType):
        supported_versions = versioning.AsdfSpec('>=1.2.0')

    assert TestType0.incompatible_version('1.1.0') == True
    assert TestType0.incompatible_version('1.2.0') == False
    assert TestType0.incompatible_version('2.0.1') == False

    class TestType1(types.CustomType):
        supported_versions = versioning.AsdfVersion('1.0.0')

    assert TestType1.incompatible_version('1.0.0') == False
    assert TestType1.incompatible_version('1.1.0') == True

    class TestType2(types.CustomType):
        supported_versions = '1.0.0'

    assert TestType2.incompatible_version('1.0.0') == False
    assert TestType2.incompatible_version('1.1.0') == True

    class TestType3(types.CustomType):
        # This doesn't make much sense, but it's just for the sake of example
        supported_versions = ['1.0.0', versioning.AsdfSpec('>=2.0.0')]

    assert TestType3.incompatible_version('1.0.0') == False
    assert TestType3.incompatible_version('1.1.0') == True
    assert TestType3.incompatible_version('2.0.0') == False
    assert TestType3.incompatible_version('2.0.1') == False

    class TestType4(types.CustomType):
        supported_versions = ['1.0.0', versioning.AsdfVersion('1.1.0')]

    assert TestType4.incompatible_version('1.0.0') == False
    assert TestType4.incompatible_version('1.0.1') == True
    assert TestType4.incompatible_version('1.1.0') == False
    assert TestType4.incompatible_version('1.1.1') == True

    class TestType5(types.CustomType):
        supported_versions = \
            [versioning.AsdfSpec('<1.0.0'), versioning.AsdfSpec('>=2.0.0')]

    assert TestType5.incompatible_version('0.9.9') == False
    assert TestType5.incompatible_version('2.0.0') == False
    assert TestType5.incompatible_version('2.0.1') == False
    assert TestType5.incompatible_version('1.0.0') == True
    assert TestType5.incompatible_version('1.1.0') == True

    with pytest.raises(ValueError):
        class TestType6(types.CustomType):
            supported_versions = 'blue'
    with pytest.raises(ValueError):
        class TestType7(types.CustomType):
            supported_versions = ['1.1.0', '2.2.0', 'blue']

def test_supported_versions():
    class CustomFlow:
        def __init__(self, c=None, d=None):
            self.c = c
            self.d = d

    class CustomFlowType(types.CustomType):
        version = '1.1.0'
        supported_versions = [(1,0,0), versioning.AsdfSpec('>=1.1.0')]
        name = 'custom_flow'
        organization = 'nowhere.org'
        standard = 'custom'
        types = [CustomFlow]

        @classmethod
        def from_tree(cls, tree, ctx):
            # Convert old schema to new CustomFlow type
            if cls.version == '1.0.0':
                return CustomFlow(c=tree['a'], d=tree['b'])
            else:
                return CustomFlow(**tree)

        @classmethod
        def to_tree(cls, data, ctx):
            if cls.version == '1.0.0':
                return dict(a=data.c, b=data.d)
            else:
                return dict(c=data.c, d=data.d)

    class CustomFlowExtension(CustomExtension):
        @property
        def types(self):
            return [CustomFlowType]

    new_yaml = """
flow_thing:
  !<tag:nowhere.org:custom/custom_flow-1.1.0>
    c: 100
    d: 3.14
"""
    old_yaml = """
flow_thing:
  !<tag:nowhere.org:custom/custom_flow-1.0.0>
    a: 100
    b: 3.14
"""
    new_buff = helpers.yaml_to_asdf(new_yaml)
    new_data = asdf.open(new_buff, extensions=CustomFlowExtension())
    assert type(new_data.tree['flow_thing']) == CustomFlow

    old_buff = helpers.yaml_to_asdf(old_yaml)
    old_data = asdf.open(old_buff, extensions=CustomFlowExtension())
    assert type(old_data.tree['flow_thing']) == CustomFlow

def test_unsupported_version_warning():
    class CustomFlow:
        pass

    class CustomFlowType(types.CustomType):
        version = '1.0.0'
        supported_versions = [(1,0,0)]
        name = 'custom_flow'
        organization = 'nowhere.org'
        standard = 'custom'
        types = [CustomFlow]

    class CustomFlowExtension(CustomExtension):
        @property
        def types(self):
            return [CustomFlowType]

    yaml = """
flow_thing:
  !<tag:nowhere.org:custom/custom_flow-1.1.0>
    c: 100
    d: 3.14
"""
    buff = helpers.yaml_to_asdf(yaml)

    with pytest.warns(AsdfConversionWarning, match="Version 1.1.0 of tag:nowhere.org:custom/custom_flow is not compatible"):
        asdf.open(buff, extensions=CustomFlowExtension())

def test_tag_without_schema(tmpdir):

    tmpfile = str(tmpdir.join('foo.asdf'))

    class FooType(types.CustomType):
        name = 'foo'

        def __init__(self, a, b):
            self.a = a
            self.b = b

        @classmethod
        def from_tree(cls, tree, ctx):
            return cls(tree['a'], tree['b'])

        @classmethod
        def to_tree(cls, node, ctx):
            return dict(a=node.a, b=node.b)

        def __eq__(self, other):
            return self.a == other.a and self.b == other.b

    class FooExtension:
        @property
        def types(self):
            return [FooType]

        @property
        def tag_mapping(self):
            return []

        @property
        def url_mapping(self):
            return []

    foo = FooType('hello', 42)
    tree = dict(foo=foo)

    with pytest.warns(AsdfWarning, match="Unable to locate schema file"):
        with asdf.AsdfFile(tree, extensions=FooExtension()) as af:
            af.write_to(tmpfile)

    with pytest.warns(AsdfWarning, match="Unable to locate schema file"):
        with asdf.AsdfFile(tree, extensions=FooExtension()) as ff:
            assert isinstance(ff.tree['foo'], FooType)
            assert ff.tree['foo'] == tree['foo']

def test_custom_reference_cycle(tmpdir):
    f1 = FractionWithInverse(3, 5)
    f2 = FractionWithInverse(5, 3)
    f1.inverse = f2
    f2.inverse = f1
    tree = {"fraction": f1}

    path = str(tmpdir.join("with_inverse.asdf"))

    with asdf.AsdfFile(tree, extensions=FractionWithInverseExtension()) as af:
        af.write_to(path)

    with asdf.open(path, extensions=FractionWithInverseExtension()) as af:
        assert af["fraction"].inverse.inverse is af["fraction"]

max_stars_repo_name: zanecodes/asdf

max_issues_repo_name: zanecodes/asdf

max_forks_repo_name: zanecodes/asdf

reason: from numpy


text:
#!/usr/bin/env python
"""
.. module:: test_all
   :synopsis: Nose tests of the lcad language.

.. moduleauthor:: Hazen Babcock

"""

import math
import nose
import numbers
import numpy

import opensdraw.lcad_language.belt as belt
import opensdraw.lcad_language.chain as chain
import opensdraw.lcad_language.curve as curve
import opensdraw.lcad_language.interpreter as interpreter
import opensdraw.lcad_language.lcadExceptions as lcadExceptions
import opensdraw.lcad_language.lexerParser as lexerParser
import opensdraw.lcad_language.lcadTypes as lcadTypes
import opensdraw.lcad_language.pulleySystem as pulleySystem

def exe(string):
    """
    Wrap interpreter call for convenience.
    """
    lenv = interpreter.LEnv(add_built_ins = True)
    model = interpreter.Model()
    ast = lexerParser.parse(string, "test")
    interpreter.createLexicalEnv(lenv, ast)
    sym = interpreter.interpret(model, ast)
    return interpreter.getv(sym)

## Symbols

# t (true)
def test_t():
    assert exe("(if t 0 1)") == 0

# nil (false)
def test_nil():
    assert exe("(if nil 0 1)") == 1

# e
def test_e():
    assert exe("e") == math.e

# pi
def test_pi():
    assert exe("pi") == math.pi

# time-index
def text_time_index():
    assert exe("time-index") == 0

## Comparison Functions.

# equal
def test_eq_1():
    assert exe("(if (= 1 1) 0 1)") == 0

def test_eq_2():
    assert exe("(if (= 1 0) 0 1)") == 1

# gt
def test_gt_1():
    assert exe("(if (> 1 2) 0 1)") == 1

def test_gt_2():
    assert exe("(if (> 1 0) 0 1)") == 0

# lt
def test_lt_1():
    assert exe("(if (< 1 2) 0 1)") == 0

def test_lt_2():
    assert exe("(if (< 1 0) 0 1)") == 1

# le
def test_le_1():
    assert exe("(if (<= 1 1) 0 1)") == 0

def test_le_2():
    assert exe("(if (<= 1 0) 0 1)") == 1

# ge
def test_ge_1():
    assert exe("(if (>= 1 2) 0 1)") == 1

def test_ge_2():
    assert exe("(if (>= 1 1) 0 1)") == 0

# ne
def test_ne_1():
    assert exe("(if (!= 1 2) 0 1)") == 0

def test_ne_2():
    assert exe("(if (!= 1 1) 0 1)") == 1

## coreFunctions

# append
def test_append_1():
    assert exe("(def m (list 1)) (append m 2) (aref m 1)") == 2

def test_append_2():
    assert exe("(def m (list 1)) (append m 2 3) (aref m 2)") == 3

# aref
def test_aref_1():
    assert exe("(aref (list 1 2 3) 1)") == 2

def test_aref_2():
    assert exe("(def x (list 1 2 3)) (set (aref x 1) 4) (aref x 1)") == 4

def test_aref_3():
    assert exe("(def x (vector 1 2 3)) (set (aref x 1) 4) (aref x 1)") == 4

def test_aref_4():
    assert exe("(def m (matrix (list 1 2 3 4 5 6 7 8 9 10 11 12))) (set (aref m 1 1) 0) (aref m 1 1)") == 0

@nose.tools.raises(lcadExceptions.LCadException)
def test_aref_5():
    exe("(aref (list 0 0 0 0 0 0) 1 1)")

@nose.tools.raises(lcadExceptions.OutOfRangeException)
def test_aref_6():
    exe("(aref (list 1 2 3) 5)")

@nose.tools.raises(lcadExceptions.LCadException)
def test_aref_7():
    exe("(aref (vector 0 0 0) 1 2)")

@nose.tools.raises(lcadExceptions.LCadException)
def test_aref_8():
    exe("(aref (matrix (list 0 0 0 0 0 0)) 5 5)")

@nose.tools.raises(lcadExceptions.WrongTypeException)
def test_aref_9():
    exe("(aref 1 5 5)")

# block
def test_block_1():
    assert exe("(def fn (block (def x 5) (def inc-x () (+ x 1)) inc-x)) (fn)") == 6

def test_block_2():
    assert exe("(def fn () (block)) (if (fn) 0 1)") == 1

# concatenate
def test_concatenate_1():
    assert exe("(concatenate \"as\" 1)") == "as1"

# cond
def test_cond_1():
    assert exe("(def x 2) (cond ((= x 1) 2) ((= x 2) 3) ((= x 3) 4) (t 5))") == 3

# copy
def test_copy_1():
    assert exe("(def a (list 1 2 3) b (copy a)) (set (aref b 0) 5) (aref a 0)") == 1

# def
def test_def_1():
    assert exe("(def x 15) x") == 15

def test_def_2():
    assert exe("(def x 15 y 20) y") == 20

def test_def_3():
    assert exe("(def incf (x) (+ x 1)) (incf 2)") == 3

def test_def_4():
    assert exe("(def incf (x :y 0) (+ x y 1)) (incf 1)") == 2

def test_def_5():
    assert exe("(def incf (x :y 0) (+ x y 1)) (incf 1 :y 2)") == 4

def test_def_6():
    assert exe("(def incf (:y 0) (+ y 1)) (incf :y 2)") == 3

def test_def_7():
    assert exe("(def incf (:y 0) (+ y 1)) (incf)") == 1

@nose.tools.raises(lcadExceptions.NumberArgumentsException)
def test_def_8():
    exe("(def a)")

@nose.tools.raises(lcadExceptions.NumberArgumentsException)
def test_def_9():
    exe("(def a 1 b 2 c)")

@nose.tools.raises(lcadExceptions.CannotSetException)
def test_def_10():
    exe("(def 1 1)")

@nose.tools.raises(lcadExceptions.CannotOverrideBuiltIn)
def test_def_11():
    exe("(def t 1)")

@nose.tools.raises(lcadExceptions.SymbolAlreadyExists)
def test_def_11():
    exe("(def x 1) (def x 2)")

@nose.tools.raises(lcadExceptions.NotAFunctionException)
def test_def_12():
    exe("(def a 1) (a)")

@nose.tools.raises(TypeError)
def test_def_13():
    exe("(def fn (x :y 2) (+ x y)) (fn 1 2)")

@nose.tools.raises(lcadExceptions.KeywordException)
def test_def_14():
    exe("(def fn (x :y 2) (+ x y)) (fn 1 y 2)")

@nose.tools.raises(lcadExceptions.UnknownKeywordException)
def test_def_15():
    exe("(def fn (x :y 2) (+ x y)) (fn 1 :z 3)")

# for
def test_for_1():
    assert exe("(def x 0) (for (i 10) (set x (+ 1 x))) x") == 10

def test_for_2():
    assert exe("(def x 0) (for (i 1 11) (set x (+ 1 x))) x") == 10

def test_for_3():
    assert exe("(def x 0) (for (i 1 0.1 2) (set x (+ 1 x))) x") == 10

def test_for_4():
    assert exe("(def x 0) (for (i (list 1 2 3)) (set x (+ i x))) x") == 6

@nose.tools.raises(lcadExceptions.NumberArgumentsException)
def test_for_5():
    exe("(for 10)")

@nose.tools.raises(lcadExceptions.LCadException)
def test_for_6():
    exe("(for i 1)")

@nose.tools.raises(lcadExceptions.NumberArgumentsException)
def test_for_7():
    exe("(for (i) 1)")

@nose.tools.raises(lcadExceptions.NumberArgumentsException)
def test_for_8():
    exe("(for (i 1 10 20 30) 1)")

@nose.tools.raises(lcadExceptions.LCadException)
def test_for_9():
    exe("(for (1 10) 1)")

# if
def test_if_1():
    assert exe("(if t 1 2)") == 1

def test_if_2():
    assert exe("(if (= 1 2) 1 2)") == 2

def test_if_3():
    assert exe("(if (if (= 1 2) t nil) 1 2)") == 2

def test_if_4():
    assert exe("(if (if (= 1 2) t) 1 2)") == 2

@nose.tools.raises(lcadExceptions.NumberArgumentsException)
def test_if_5():
    exe("(if t)")

@nose.tools.raises(lcadExceptions.NumberArgumentsException)
def test_if_6():
    exe("(if t 0 1 2)")

# import
def test_import_1():
    assert exe("(import mod) (mod:fn)") == math.pi

def test_import_2():
    assert exe("(import mod :local) (fn)") ==  math.pi

@nose.tools.raises(lcadExceptions.NumberArgumentsException)
def test_import_3():
    exe("(import)")

# is-main
def test_is_main_1():
    assert exe("(if (is-main) 1 2)") == 1

# lambda
def test_lambda_1():
    assert exe("(def fn (lambda (x) (+ x 1))) (fn 1)") == 2

@nose.tools.raises(lcadExceptions.NumberArgumentsException)
def test_lambda_2():
    exe("(lambda (x))")

# len
def test_len_1():
    assert exe("(len (list 1 2 3))") == 3

# list
def test_list_1():
    assert exe("(def x (list 1 2 3)) (aref x 0)") == 1

def test_list_2():
    assert exe("(list 1 2 3)")[1] == 2

# print
def test_print_1():
    assert exe("(print \"123\")") == "123"

# pyimport
def test_pyimport_1():
    assert exe("(pyimport pyimp1) (plus 1 1)") == 2

@nose.tools.raises(lcadExceptions.NumberArgumentsException)
def test_pyimport_2():
    exe("(pyimport)")

# set
def test_set_1():
    assert exe("(def x 10) (set x 15)") == 15

def test_set_2():
    assert exe("(def x 10 y 10) (set x 15 y 20)") == 20

def test_set_3():
    assert exe("(def fn () 1) (def x 2) (set x fn) (x)") == 1

@nose.tools.raises(lcadExceptions.NumberArgumentsException)
def test_set_4():
    exe("(set x)")

@nose.tools.raises(lcadExceptions.CannotSetException)
def test_set_5():
    exe("(set 1 2)")

@nose.tools.raises(lcadExceptions.CannotOverrideBuiltIn)
def test_set_6():
    exe("(set t nil)")

# while
def test_while_1():
    assert exe("(def x 0) (while (< x 9) (set x (+ 2 x))) x") == 10

@nose.tools.raises(lcadExceptions.NumberArgumentsException)
def test_while_2():
    exe("(while t)")

## Geometry Functions.

# cross product
def test_cross_product_1():
    assert exe("(aref (cross-product (vector 1 0 0) (vector 0 1 0)) 2)") == 1

def test_cross_product_2():
    assert exe("(aref (cross-product (vector 1 0 0) (vector 0 2 0) nil) 2)") == 2

# dot product
def test_dot_product_1():
    assert exe("(dot-product (vector 1 0 0) (vector 1 0 0))") == 1

def test_dot_product_2():
    assert exe("(dot-product (vector 1 0 0) (vector 0 1 0))") == 0

# matrix
def test_matrix_1():
    assert isinstance(exe("(matrix (list 1 2 3 1 2 3 1 2 3 1 2 3))"), lcadTypes.LCadMatrix)

def test_matrix_2():
    assert isinstance(exe("(matrix (list 1 2 3 1 2 3))"), lcadTypes.LCadMatrix)

def test_matrix_3():
    assert isinstance(exe("(matrix (list (list 1 2 3) (list 1 2 3) (list 1 2 3) (list 1 2 3)))"), lcadTypes.LCadMatrix)

def test_matrix_4():
    assert isinstance(exe("(matrix (matrix (list 1 2 3 1 2 3 1 2 3 1 2 3)))"), lcadTypes.LCadMatrix)

# mirror
def test_mirror_1():
    assert exe("(mirror (list 1 (if t 0 1) (if nil 0 1)) 1)") == 1

def test_mirror_2():
    assert exe("(mirror (vector 1 0 0) 1)") == 1

# rotate
def test_rotate_1():
    assert exe("(rotate (list 1 2 3) 1)") == 1

def test_rotate_2():
    assert exe("(rotate (vector 1 2 3) 1)") == 1

# scale
def test_scale_1():
    assert exe("(scale (list 1 2 3) 1)") == 1

def test_scale_2():
    assert exe("(scale (vector 1 2 3) 1)") == 1

# transform
def test_transform_1():
    assert exe("(transform (list 1 2 3 1 2 3 1 2 3 1 2 3) 1)") == 1

def test_transform_2():
    assert exe("(transform (list 1 2 3 1 2 3) 1)") == 1

def test_transform_3():
    assert exe("(transform (matrix (list 1 2 3 1 2 3 1 2 3 1 2 3)) 1)") == 1

# translate
def test_translate_1():
    assert exe("(translate (list 1 2 3) 1)") == 1

def test_translate_1():
    assert exe("(translate (vector 1 2 3) 1)") == 1

# vector
def test_vector_1():
    assert isinstance(exe("(vector 1 2 3)"), lcadTypes.LCadVector)

@nose.tools.raises(lcadExceptions.NumberArgumentsException)
def test_vector_2():
    exe("(vector 1)")

@nose.tools.raises(lcadExceptions.WrongTypeException)
def test_vector_3():
    exe("(vector 1 2 'a')")

@nose.tools.raises(lcadExceptions.WrongTypeException)
def test_vector_4():
    exe("(vector 1 2 3 4 'a')")

## Logic Functions.

# and
def test_and_1():
    assert exe("(if (and (< 1 2) (< 2 3)) 0 1)") == 0

def test_and_2():
    assert exe("(if (and (< 1 2) (> 2 3)) 0 1)") == 1

# or
def test_or_1():
    assert exe("(if (or (> 1 2) (< 2 3)) 0 1)") == 0

def test_or_2():
    assert exe("(if (or (> 1 2) (> 2 3)) 0 1)") == 1

# not
def test_not_1():
    assert exe("(if (not t) 0 1)") == 1

def test_not_2():
    assert exe("(if (not nil) 0 1)") == 0

def test_not_3():
    assert exe("(if (not ()) 0 1)") == 0

## Math Functions.

# basic math
def test_math_1():
    assert exe("(+ 1 1)") == 2

def test_math_2():
    assert exe("(- 1 1)") == 0

def test_math_3():
    assert exe("(* 2 2)") == 4

def test_math_4():
    assert exe("(/ 4 2)") == 2

def test_math_5():
    assert exe("(% 11 2)") == 1

def test_math_6():
    assert exe("(* (matrix (list 1 1 1 3 0 0 0 3 0 0 0 3)) (vector 1 2 3 1))")[1] == 7

def test_math_7():
    assert isinstance(exe("(* (matrix (list 1 1 1 3 0 0 0 3 0 0 0 3)) (vector 1 2 3 1))"), lcadTypes.LCadVector)

def test_math_8():
    assert isinstance(exe("(* (matrix (list 1 2 3 1 2 3)) (matrix (list 1 2 3 1 2 3)))"), lcadTypes.LCadMatrix)

def test_math_9():
    assert exe("(* (vector 1 2 3) (vector 1 2 3))")[1] == 4

def test_math_10():
    assert exe("(/ (vector 1 2 3) (vector 1 2 3))")[1] == 1

def test_math_11():
    assert exe("(+ (vector 1 2 3) (vector 1 2 3))")[1] == 4

def test_math_12():
    assert exe("(- (vector 1 2 3) (vector 1 2 3))")[1] == 0

def test_math_13():
    assert exe("(- 1)") == -1

def test_math_14():
    assert exe("(abs 2)") == 2

def test_math_15():
    assert exe("(abs -2)") == 2

# python math module
def test_py_math_1():
    assert int(round(exe("(cos 0)"))) == 1

def test_py_math_2():
    assert int(round(exe("(cos (/ pi 2))"))) == 0

def test_py_math_3():
    assert int(round(exe("(sin 0)"))) == 0

def test_py_math_4():
    assert int(round(exe("(sin (/ pi 2))"))) == 1

## Part Functions.

# group
def test_group_1():
    assert exe("(group \"adsf\" 1)") == 1

# header
def test_header_1():
    assert exe("(header \"asdf\")") == "asdf"

# line
def test_line_1():
    assert exe("(line (list 1 2 3) (list 4 5 6)) 1") == 1

def test_line_2():
    assert exe("(line (vector 1 2 3) (vector 4 5 6)) 1") == 1

# optional line
def test_optional_line_1():
    assert exe("(optional-line (list 1 2 3) (list 4 5 6) (list 1 2 3) (list 4 5 6)) 1") == 1

def test_optional_line_1():
    assert exe("(optional-line (vector 1 2 3) (vector 4 5 6) (vector 1 2 3) (vector 4 5 6)) 1") == 1

# part
def test_part_1():
    assert exe("(part '1234' 5) 1") == 1

# quadrilateral
def test_quadrilateral_1():
    assert exe("(quadrilateral (list 1 2 3) (list 4 5 6) (list 1 2 3) (list 4 5 6)) 1") == 1

def test_quadrilateral_2():
    assert exe("(quadrilateral (vector 1 2 3) (vector 4 5 6) (vector 1 2 3) (vector 4 5 6)) 1") == 1

# triangle
def test_triangle_1():
    assert exe("(triangle (list 1 2 3) (list 4 5 6) (list 7 8 9)) 1") == 1

def test_triangle_1():
    assert exe("(triangle (vector 1 2 3) (vector 4 5 6) (vector 7 8 9)) 1") == 1

# Random Number Functions.
def test_rand_seed_1():
    assert isinstance(exe("(rand-seed 10)"), numbers.Number)

def test_rand_choice_1():
    assert isinstance(exe("(rand-choice (list 1 2 3))"), numbers.Number)

def test_rand_choice_2():
    assert isinstance(exe("(def a () 10) (def b () 20) (def c (rand-choice (list a b))) (c)"), numbers.Number)

def test_rand_gauss_1():
    assert isinstance(exe("(rand-gauss)"), numbers.Number)

def test_rand_gauss_2():
    assert isinstance(exe("(rand-gauss 1 2)"), numbers.Number)

def test_rand_integer_1():
    assert isinstance(exe("(rand-integer 0 10)"), numbers.Number)

def test_rand_uniform_1():
    assert isinstance(exe("(rand-uniform)"), numbers.Number)

def test_rand_uniform_2():
    assert isinstance(exe("(rand-uniform 1 10)"), numbers.Number)

## Type Functions.
def test_is_boolean_1():
    assert exe("(if (boolean? nil) 1 0)") == 1

def test_is_boolean_2():
    assert exe("(if (boolean? 1) 1 0)") == 0

def test_is_matrix_1():
    assert exe("(if (matrix? (matrix (list 0 0 0 0 0 0))) 1 0)") == 1

def test_is_matrix_2():
    assert exe("(if (matrix? (vector 0 0 0)) 1 0)") == 0

def test_is_number_1():
    assert exe("(if (number? 1) 1 0)") == 1

def test_is_number_2():
    assert exe("(if (number? 'a') 1 0)") == 0

def test_is_string_1():
    assert exe("(if (string? 1) 1 0)") == 0

def test_is_string_2():
    assert exe("(if (string? 'a') 1 0)") == 1

def test_is_vector_1():
    assert exe("(if (vector? (matrix (list 0 0 0 0 0 0))) 1 0)") == 0

def test_is_vector_2():
    assert exe("(if (vector? (vector 0 0 0)) 1 0)") == 1

## Library Functions.

# belt
def test_belt_1():
    assert exe("(belt (list (list (list 0 0 0) (list 0 0 1) 1.0 1) (list (list 4 0 0) (list 0 0 1) 1.5 1))) 1") == 1

def test_belt_2():
    assert exe("(belt (list (list (list 0 0 0) (list 0 0 1) 1.0 1) (list (list 4 0 0) (list 0 0 1) 1.5 1)) :continuous nil) 1") == 1

@nose.tools.raises(belt.NumberPulleysException)
def test_belt_3():
    exe("(belt (list (list (list 0 0 0) (list 0 0 1) 1.0 1)))")

@nose.tools.raises(lcadExceptions.WrongTypeException)
def test_belt_4():
    exe("(belt (list (vector 1 1 1 1) (vector 1 1 1 1)))")

@nose.tools.raises(belt.PulleyException)
def test_belt_5():
    exe("(belt (list (list (list 0 0 0) 1.0 1) (list (list 4 0 0) (list 0 0 1) 1.5 1)))")

@nose.tools.raises(lcadExceptions.LCadException)
def test_belt_6():
    exe("(belt (list (list (list 0 0) (list 0 0 1) 1.0 1) (list (list 4 0 0) (list 0 0 1) 1.5 1)))")

@nose.tools.raises(lcadExceptions.LCadException)
def test_belt_7():
    exe("(belt (list (list (list 0 0 0) (list 0 0) 1.0 1) (list (list 4 0 0) (list 0 0 1) 1.5 1)))")

@nose.tools.raises(lcadExceptions.WrongTypeException)
def test_belt_8():
    exe("(belt (list (list (list 0 0 0) (list 0 0 1) t 1) (list (list 4 0 0) (list 0 0 1) 1.5 1)))")

@nose.tools.raises(lcadExceptions.WrongTypeException)
def test_belt_9():
    exe("(belt (list (list (list 0 0 0) (list 0 0 1) 1.0 'a') (list (list 4 0 0) (list 0 0 1) 1.5 1)))")

# chain
def test_chain_1():
    assert exe("(chain (list (list -4 0 1 1) (list 4 0 1 1))) 1") == 1

def test_chain_2():
    assert exe("(chain (list (list -4 0 1 1) (list 4 0 1 1)) :continuous nil) 1") == 1

@nose.tools.raises(chain.NumberSprocketsException)
def test_chain_3():
    exe("(chain (list (list -4 0 1 1)))")

@nose.tools.raises(lcadExceptions.WrongTypeException)
def test_chain_4():
    exe("(chain (list (list -4 0 1 1) (vector 4 0 1 1)))")

@nose.tools.raises(chain.SprocketException)
def test_chain_5():
    exe("(chain (list (list -4 0 1 1) (list 0 1 1)))")

@nose.tools.raises(lcadExceptions.WrongTypeException)
def test_chain_6():
    exe("(chain (list (list -4 0 1 1) (list t 0 1 1)))")

# curve
def test_curve_1():
    assert exe("(curve (list (list (list 0 0 0) (list 1 1 0) (list 0 0 1)) (list (list 5 0 0) (list 1 0 0)))) 1") == 1

def test_curve_2():
    assert exe("(def my-curve (curve (list (list (list 0 0 0) (list 1 0 0) (list 0 0 1)) (list (list 1 0 0) (list 1 0 0))))) (my-curve t)") == 1.0

def test_curve_3():
    assert exe("(curve (list (list (list 0 0 0) (list 1 1 0) (list 0 0 1)) (list (list 5 0 0) (list 1 0 0))) :auto-scale nil) 1") == 1

@nose.tools.raises(curve.NumberControlPointsException)
def test_curve_4():
    exe("(curve (list (list (list 0 0 0) (list 1 1 0) (list 0 0 1))))")

@nose.tools.raises(lcadExceptions.WrongTypeException)
def test_curve_5():
    exe("(curve (list (list (list 0 0 0) (list 1 1 0) (list 0 0 1)) (vector 0 1 2)))")

@nose.tools.raises(curve.ControlPointException)
def test_curve_6():
    exe("(curve (list (list (list 0 0 0) (list 0 0 1)) (list (list 5 0 0) (list 1 0 0))))")

@nose.tools.raises(curve.ControlPointException)
def test_curve_7():
    exe("(curve (list (list (list 0 0 0) (list 1 1 0) (list 0 0 1)) (list (list 5 0 0))))")

@nose.tools.raises(curve.ControlPointException)
def test_curve_8():
    exe("(curve (list (list (list 0 0) (list 1 1 0) (list 0 0 1)) (list (list 5 0 0) (list 1 0 0))))")

@nose.tools.raises(lcadExceptions.WrongTypeException)
def test_curve_9():
    exe("(curve (list (list (list t 0 0) (list 1 1 0) (list 0 0 1)) (list (list 5 0 0) (list 1 0 0))))")

@nose.tools.raises(TypeError)
def test_curve_10():
    exe("(curve (list (list t (list 1 1 0) (list 0 0 1)) (list (list 5 0 0) (list 1 0 0))))")

@nose.tools.raises(curve.TangentException)
def test_curve_11():
    exe("(curve (list (list (list 0 0 0) (list 0 0 0) (list 0 0 1)) (list (list 5 0 0) (list 1 0 0))))")

# pulley system.
def test_pulley_system_1():
    assert exe("(pulley-system (list (list (list 0 0 0) (list 0 0 1) 5.0 1 5.0 1 50) (list (list 20 0 1) (list 0 0 1) 5.0 1) (list (list -1 0 0) \"tangent\"))) 1") == 1

def test_pulley_system_2():
    assert exe("(pulley-system (list (list (list 0 0 0) (list 0 0 1) 5.0 1 5.0 1 50) (list (list 20 0 0) \"point\"))) 1") == 1

@nose.tools.raises(pulleySystem.NumberPulleysException)
def test_pulley_system_3():
    exe("(pulley-system (list (list (list 0 0 0) (list 0 0 1) 5.0 1 5.0 1 50)))")

@nose.tools.raises(lcadExceptions.WrongTypeException)
def test_pulley_system_4():
    exe("(pulley-system (list (vector 1 2 3) (list (list 20 0 0) \"point\")))")

@nose.tools.raises(pulleySystem.DrumException)
def test_pulley_system_5():
    exe("(pulley-system (list (list (list 0 0 0) (list 0 0 1) 5.0 1 5.0 1) (list (list 20 0 0) \"point\")))")

@nose.tools.raises(lcadExceptions.WrongTypeException)
def test_pulley_system_6():
    exe("(pulley-system (list (list (list 0 0 0) (list 0 0 1) t 1 5.0 1 50) (list (list 20 0 0) \"point\")))")

@nose.tools.raises(lcadExceptions.WrongTypeException)
def test_pulley_system_7():
    exe("(pulley-system (list (list (list 0 0 0) (list 0 0 1) 5.0 1 5.0 1 50) (vector 0 0 0)))")

@nose.tools.raises(pulleySystem.EndPointException)
def test_pulley_system_8():
    exe("(pulley-system (list (list (list 0 0 0) (list 0 0 1) 5.0 1 5.0 1 50) (list (list 20 0 0))))")

@nose.tools.raises(lcadExceptions.WrongTypeException)
def test_pulley_system_9():
    exe("(pulley-system (list (list (list 0 0 0) (list 0 0 1) 5.0 1 5.0 1 50) (list (list 20 0 0) t)))")

@nose.tools.raises(pulleySystem.EndPointTypeException)
def test_pulley_system_10():
    exe("(pulley-system (list (list (list 0 0 0) (list 0 0 1) 5.0 1 5.0 1 50) (list (list 20 0 0) \"foo\")))")

# spring
def test_spring_1():
    assert exe("(spring 40 10 1 10) 1") == 1

def test_spring_2():
    assert exe("(spring 40 10 1 10 1) 1") == 1

max_stars_repo_name: HazenBabcock/openldraw

max_issues_repo_name: HazenBabcock/openldraw

max_forks_repo_name: HazenBabcock/openldraw


text:
import os
import torch.utils.data as data
#using the rotating building data...
#bu_filename = 'test_5010r'
import geopandas as gpd
from shapely import affinity
from shapely.geometry import Polygon
import numpy as np
import copy
from new_pre import *
bu_filename = 'test_5010'     #'test_5010r'
hparams = {
    # preprocessing params
    # 曲线变直线方法参数
    'if_simplify': True
    , 'tor_dist': 0.1
    , 'tor_cos': 0.99
    , 'if_scale_y': False

    , 'scale_type': 'final'
    , 'seq_length': 32  #
    , 'rotate_type': 'equal'  #
    , 'rotate_length': 1  #
    , 'norm_type': 'minmax'  #

    # model params
    , 'GPU': True  #
    , 'epochs': 400  #
    , 'optimizer': 'rmsprop'  # 均方根传播法
    , 'z_size': 128  #
    , 'rnn_type': 'lstm'  # lstm or gru
}
bu_shape = gpd.read_file('./data/' + bu_filename + '.shp', encode='utf-8')
bu_use = copy.deepcopy(bu_shape)  # deepcopy的结果不会随原数据的改变而改变
bu_mbr, bu_use = get_shape_mbr(bu_use)
bu_use = get_shape_normalize_final(bu_use, hparams['if_scale_y'])  # hparams在55行自定义参数 False
if hparams['if_simplify']:  # True
    bu_use = get_shape_simplify(bu_use, hparams['tor_dist'], hparams['tor_cos'], simplify_type=0)  # 0.1 0.99
bu_use = reset_start_point(bu_use)  # 重设多边形起始点
bu_node = get_node_features(bu_use)
bu_line = get_line_features_final(bu_node, hparams['seq_length'])
bu_detail = get_inter_features(bu_line)  # 得到插入特征
# outputpath = "./data/bu_detail.csv"
# bu_detail.to_csv(outputpath, sep=',', index=True, header=True)
bu_detail = get_neat_features(bu_detail, hparams['seq_length'], hparams['rotate_length'])  # 处理，得到平滑特征
# outputpath = "./data/bu_detail_neat.csv"
# bu_detail.to_csv(outputpath, sep=',', index=True, header=True)
test_x = np.array(bu_detail['xs']).reshape(5010*hparams['seq_length'], 1)
test_y = np.array(bu_detail['ys']).reshape(5010*hparams['seq_length'], 1)
test_x_y = np.concatenate((test_x, test_y), axis=1)
test_xy_reshape = test_x_y.reshape(5010, hparams['seq_length'], 2)

class ModelNet40DataSet(data.Dataset):
    def __init__(self, train=True):
        label = np.genfromtxt("data/label.content", dtype=np.int64).reshape(-1,1)
        modelnet_data = np.zeros([0,hparams['seq_length'],2], dtype=np.float64)
        modelnet_label = np.zeros([0,1], np.float64)
        modelnet_data = np.concatenate([modelnet_data, test_xy_reshape],axis=0)
        modelnet_label = np.concatenate([modelnet_label,label],axis=0)
        self.point_cloud = modelnet_data
        self.label = modelnet_label
    def __getitem__(self, item):
        return self.point_cloud[item], self.label[item]

    def __len__(self):
        return self.label.shape[0]

if __name__ == '__main__':
    modelnet_train = ModelNet40DataSet(train=True)
    modelnet_test = ModelNet40DataSet(train=False)

    data, label = modelnet_train[0]
    print('type: {}'.format(type(data)))
    print('type: {}'.format(type(label)))
    print(label)


text:
# -*- coding: utf-8 -*-
"""
Created on Sat Aug 24 14:44:07 2013

@author: steve
"""

from nose.tools import assert_true
import numpy as np
import scipy as sp

import hiive.mdptoolbox.example
from hiive import mdptoolbox

class BaseTestMDP(object):
    small_P, small_R = mdptoolbox.example.small()

def assert_sequence_almost_equal(a, b, spacing=10e-12):
    assert_true(all(abs(a[k] - b[k]) < spacing for k in range(len(a))))

STATES = 10
ACTIONS = 3
SMALLNUM = 10e-12

# np.arrays
P_small, R_small = mdptoolbox.example.small()

P_sparse = np.empty(2, dtype=object)
P_sparse[0] = sp.sparse.csr_matrix(P_small[0])
P_sparse[1] = sp.sparse.csr_matrix(P_small[1])

P_forest, R_forest = mdptoolbox.example.forest()

P_forest_sparse, R_forest_sparse = mdptoolbox.example.forest(S=STATES,
                                                             is_sparse=True)

np.random.seed(0)
P_rand, R_rand = mdptoolbox.example.rand(STATES, ACTIONS)

np.random.seed(0)
P_rand_sparse, R_rand_sparse = mdptoolbox.example.rand(STATES, ACTIONS,
                                                       is_sparse=True)

max_stars_repo_name: rmclimber/hiivemdptoolbox

max_issues_repo_name: rmclimber/hiivemdptoolbox

max_forks_repo_name: rmclimber/hiivemdptoolbox

reason: import numpy,import scipy