I am very new to python coding and have been struggling with this for a while I cannot see the reason for this error! Any help will be greatly appreciated.

Here is the code that is causing this error:

import arcpy
import numpy as np
import math
from arcpy.sa import *

# Get the input arguements

# Input DEM
input_raster = arcpy.GetParameterAsText(0)

# Output 
output_raster = arcpy.GetParameterAsText(1)

# Input Discharge
mean_discharge = float(arcpy.GetParameterAsText(2))

# Input number of timesteps
input_timestep = int(arcpy.GetParameterAsText(3))

#Other stuff to put outside while loop

# The below text takes the input raster and calculates the bottom left corner
extent_xmin_result = arcpy.GetRasterProperties_management(input_raster, "LEFT")
extent_xmin = float(extent_xmin_result.getOutput(0))
extent_ymin_result = arcpy.GetRasterProperties_management(input_raster, "BOTTOM")
extent_ymin = float(extent_ymin_result.getOutput(0))

# Turns the corner into a point
bottom_left_corner = arcpy.Point(extent_xmin, extent_ymin)

arcpy.overwriteOutput = True

arcpy.env.workspace = r"C:\SedMod\SedMod\Workspace.gdb"
arcpy.env.scratchWorkspace = r"C:\SedMod\SedMod\Workspace.gdb"

actual_timestep = 0

while actual_timestep <= input_timestep: # Slope calculation

# Set local variables
out_measurement = "PERCENT_RISE" 

# Check out the ArcGIS Spatial Analyst extension license
arcpy.CheckOutExtension("Spatial")

#Execute Slope
out_slope = Slope(input_raster, out_measurement)

# Convert to gradient
slope_percentage = out_slope / 100

#Save slope
save_location = r"C:/SedMod/SedMod/Workspace.gdb/test" + str(actual_timestep)
slope_percentage.save(save_location)

#Execute Flow Direction

output_flowdirection = FlowDirection(input_raster)

#Execute Drainage Area

# Method for describing the raster and then using the cellsize to determine the Drainage factor for calculating Drainage area
desc_raster = arcpy.Describe(input_raster)
raster_cell_width = desc_raster.MeanCellWidth # Cell width
raster_cell_height = desc_raster.MeanCellHeight # Cell height
spatial_reference = desc_raster.spatialReference # Spatial reference
multi_factor = (raster_cell_height / 1000) * (raster_cell_width / 1000)
            
# Run the flow direction
output_flowaccumulation = FlowAccumulation(output_flowdirection)
output_drainage_area = output_flowaccumulation * multi_factor

#Execute Discharge Calculation

mean_discharge_float = float(mean_discharge)

#Calculating the discharge based on the average taken from the Avon catchment to bath (1605km2) and a mean flow of 20.762 m3/s giving a discharge per km2 of 0.01294 m3/s
discharge_raster = output_drainage_area * mean_discharge_float

#Execute StreamPower

# Set variables for calculating stream power    
density_of_water = 1000
acceleration_gravity = 9.81
    
# Stream power calculation = stream power, g is the density of water (1000 kg/m3), g is acceleration due to gravity (9.8 m/s2), discharge (m3/s) and channel slope.
unit_width_stream_power = (density_of_water * acceleration_gravity * discharge_raster * slope_percentage) / raster_cell_width

#Execute calculating Sediment Transport        
        
# Set variables for calculating dimensionless stream power    
pw = 1000 # Density of water
ps = 2650 # Density of sediment
g = 9.81 # Acceleration due to gravity 
Di = 0.01 # Sediment size
                                                      
# Dimensionless unit width stream power
dimensionless_stream_power = unit_width_stream_power / ((16186.5) * math.sqrt((1.65 * 9.81 * (Di * Di *Di ))))

# Additional variables to calculating sediment transport
timestep_hour = 3600
        
# Calculates dimensionless sediment transport       
dimensionless_sediment_transport = Con(dimensionless_stream_power < 0.25, ((dimensionless_stream_power ** 6) * 100), (0.2 * (dimensionless_stream_power ** 1.5)))
                   
# Calculates bed transport rate
bed_material_transport = dimensionless_sediment_transport * ((16186.5) * math.sqrt((1.65 * 9.81 * (Di * Di *Di ))))            
       
# Calculates it for the cell size in kg per second        
bed_material_transport_kg = (bed_material_transport * raster_cell_width)  / 9.81

# Converts it into volume m3
bed_material_transport_volume = (bed_material_transport_kg / 2650) * timestep_hour
      


# Convert the sediment transport and the flow direction rasters into Numpy arrays
sediment_transport_np = arcpy.RasterToNumPyArray(bed_material_transport_volume, '#', '#', '#', -9999)
flow_direction_np = arcpy.RasterToNumPyArray(output_flowdirection, '#', '#', '#', -9999)

[rows,cols] = flow_direction_np.shape
elevation_gain = np.zeros((rows,cols), np.float)
elevation_gain_32 = np.zeros((rows,cols), np.float)
elevation_gain_64 = np.zeros((rows,cols), np.float)
elevation_gain_128 = np.zeros((rows,cols), np.float)
elevation_gain_16 = np.zeros((rows,cols), np.float)
elevation_gain_1 = np.zeros((rows,cols), np.float)
elevation_gain_2 = np.zeros((rows,cols), np.float)
elevation_gain_4 = np.zeros((rows,cols), np.float)
elevation_gain_8 = np.zeros((rows,cols), np.float)
        
# Main body for calculating elevation change        

# Attempt 4

# Move the sediment transport in the direction for 32 direction        
for [i, j], flow in np.ndenumerate(flow_direction_np):
    try:
        if flow == 32:
            elevation_gain_32[i - 1, j - 1]  = sediment_transport_np[i, j]
    except IndexError:
            elevation_gain_32[i - 1, j - 1] = 0
            
# Move the sediment transport in the direction for 64 direction
for [i, j], flow in np.ndenumerate(flow_direction_np):
    try:
        if flow == 64:
            elevation_gain_64[i - 1, j]  = sediment_transport_np[i, j]
    except IndexError:
        elevation_gain_64[i, j] = 0

# Move the sediment transport in the direction for 128 direction                  
for [i, j], flow in np.ndenumerate(flow_direction_np):
    try:
        if flow == 128:
            elevation_gain_128[i - 1, j + 1]  = sediment_transport_np[i, j]
    except IndexError:
        elevation_gain_128[i, j] = 0                

# Move the sediment transport in the direction for 16 direction
for [i, j], flow in np.ndenumerate(flow_direction_np):
    try:
        if flow == 16:
            elevation_gain_16[i, j - 1]  = sediment_transport_np[i, j]
    except IndexError:
        elevation_gain_16[i, j] = 0

# Move the sediment transport in the direction for 1 direction
for [i, j], flow in np.ndenumerate(flow_direction_np):
    try:
        if flow == 1:
            elevation_gain_1[i, j + 1]  = sediment_transport_np[i, j]
    except IndexError:
        elevation_gain_1[i, j] = 0

# Move the sediment transport in the direction for 2 direction
for [i, j], flow in np.ndenumerate(flow_direction_np):
    try:
        if flow == 2:
            elevation_gain_2[i + 1, j + 1]  = sediment_transport_np[i, j]
    except IndexError:
        elevation_gain_2[i, j] = 0

# Move the sediment transport in the direction for 4 direction
for [i, j], flow in np.ndenumerate(flow_direction_np):
    try:
        if flow == 4:
            elevation_gain_4[i + 1, j]  = sediment_transport_np[i, j]
    except IndexError:
        elevation_gain_4[i, j] = 0

# Move the sediment transport in the direction for 8 direction
for [i, j], flow in np.ndenumerate(flow_direction_np):
    try:
        if flow == 8:
            elevation_gain_8[i + 1, j - 1]  = sediment_transport_np[i, j]
    except IndexError:
        elevation_gain_8[i, j - 1] = 0                                   

# Covert numpy arrays back into ArcGIS raster format as I do not know how to successfully add together multidimensional numpy arrays                  
elevation_gain_32_ras = arcpy.NumPyArrayToRaster(elevation_gain_32, bottom_left_corner, raster_cell_width, raster_cell_height, -9999)
elevation_gain_64_ras = arcpy.NumPyArrayToRaster(elevation_gain_64, bottom_left_corner, raster_cell_width, raster_cell_height, -9999)
elevation_gain_128_ras = arcpy.NumPyArrayToRaster(elevation_gain_128, bottom_left_corner, raster_cell_width, raster_cell_height, -9999)
elevation_gain_16_ras = arcpy.NumPyArrayToRaster(elevation_gain_16, bottom_left_corner, raster_cell_width, raster_cell_height, -9999)
elevation_gain_1_ras = arcpy.NumPyArrayToRaster(elevation_gain_1, bottom_left_corner, raster_cell_width, raster_cell_height, -9999)
elevation_gain_2_ras = arcpy.NumPyArrayToRaster(elevation_gain_2, bottom_left_corner, raster_cell_width, raster_cell_height, -9999)
elevation_gain_4_ras = arcpy.NumPyArrayToRaster(elevation_gain_4, bottom_left_corner, raster_cell_width, raster_cell_height, -9999)
elevation_gain_8_ras = arcpy.NumPyArrayToRaster(elevation_gain_8, bottom_left_corner, raster_cell_width, raster_cell_height, -9999)

# Calculate  volume in raster for whole catchment        
volume_in = elevation_gain_32_ras + elevation_gain_64_ras + elevation_gain_128_ras + elevation_gain_16_ras + elevation_gain_1_ras + elevation_gain_2_ras + elevation_gain_4_ras + elevation_gain_8_ras

#State what volume out is and the elevation           
    
net_volume = volume_in - bed_material_transport_volume

elevation_change = net_volume / (raster_cell_width * raster_cell_height)

new_elevation = input_raster + elevation_change

actual_timestep = actual_timestep + 1

new_elevation.save(output_raster)   

The error I am getting is:

Traceback (most recent call last):
  File "C:\SedMod\SedMod\Slope1.py", line 47, in <module>
    out_slope = Slope(input_raster, out_measurement, )
  File "c:\program files (x86)\arcgis\desktop10.1\arcpy\arcpy\sa\Functions.py", line 5523, in Slope
    z_factor)
  File "c:\program files (x86)\arcgis\desktop10.1\arcpy\arcpy\sa\Utils.py", line 47, in swapper
    result = wrapper(*args, **kwargs)
  File "c:\program files (x86)\arcgis\desktop10.1\arcpy\arcpy\sa\Functions.py", line 5518, in wrapper
    z_factor)
  File "c:\program files (x86)\arcgis\desktop10.1\arcpy\arcpy\geoprocessing\_base.py", line 498, in <lambda>
    return lambda *args: val(*gp_fixargs(args, True))
RuntimeError: ERROR 010240: Could not save raster dataset to C:\SedMod\SedMod\Workspace.gdb\numpy_ras6 with output format FGDBR.

Failed to execute (SedModv1).