You provide initial conditions for g[s] (equation (3)), one of them involving D(g[s]](1), but you solve an ODE which does not contain any derivative of g[s].

If g[s] must respect a differential eqution, write it and solve a 2-by-2 differential system instead of a single differential equation.

If g[s] is not governed by a differential equation, specify the expression of g[s](x) and remove all references to g[s] in the boundary conditions.
Here is an illustration of the "Known g[s] case":  Known_gs_mmcdara.mw

this proposal.mw ?

Download proposal.mw

(customize it as you want) 

See also add_on.mw

@jalal 


To get this histogram you present:

• Given the nature of your data I suppose they are already consist in a serie of elements of the type [C[k], N[k]] where
  • C[k] is the range of the kth class,
  • N[k] is the number of elements C[k] contains.

This is what I call synthetic data by opposition to raw data (see below).
 

• If it is so you have two possibilities:
  1. To build an Statistics:-Histogram you must have raw data, not synthetic data.
     So either you have these raw data or you don't.
     • If you have them, let R the list/vector of incomes, just use Statistics:-Histogram with the ad hoc customization (first plot in the attached file).
        
     • If you do not have them, you can recreate fake raw data from synthetic data so that their histogram is the one you want (second plot in the attached file).
        
  2. If you only have synthetic data and bo not necessarily want to use Statistics:-Histogram, then you can proceed as explaine in the third plot in the attached file.

Download Different_options.mw

It's likely that some typo or hidden character appeared in ode2 (bold red text).
As I never use  the document mode with 2D inputs, the attached file displays two types of fonts: so correct your own file the way I did if this annoys you.

Download PGVTN_mmcdara.mw

PS 1: using PDETools:-declare doesn't ease the task of finding where the error is and you must be doubly careful using it.
PS 2: be consistent and declare u[e] too if you consider it as a function of x (which u[e]' seems to mean, even if u[e] is later defined as a constant function).

As explained in the attached file your code does not work with my Maple 2015?
So I copy-pasted your expression of eq5 and used it as the starting point to get eq6.

Your code is in grey, mine on bold brown; the successive transformations of F (the integral term in your eq6) enableobtaoining a standard form from which it clearly appears that the orthogonality condition of Hermite polynomials (Maple uses their physical form) can be used to simplify the expression of F:

Final remark: as _C2 intervenes through its square, there are two real and opposite values of _C2 such that F(n)=1.

Download SolnforOS_dsolve_mmcdara.mw

The two derivatives give exactly the same thing a.

Nevertheless here is a way to do what you want derivatives.mw

I_would_do_that.mw

I declared prnt and clr as optional parameters (first one is boolean initialized to true, second one is string initialized to blue), I introduced a control of the color you pass, and I simplified your code a little bit.
Here are the modifications (complete code in the attached file)

spread2:=proc(
               ....
               {clr::string:="Blue"}, 
               {prnt::boolean:=true}
             )
           ....
           local col := substring(StringTools:-Capitalize(clr), 1..1);
           local COL := `if`(col="R", "Red", `if`(col="G", "Green", "Blue"))  
           ....
           if `not`(member(col, {"R", "G", "B"})) then
             error cat("allowed colors are red, green or blue, received ", clr)
           end if:

           if prnt then
             print(cat("Spread 2 [x,y] Points/Vectors wrt origin ",  COL));
           end if:
                 
           if col = "B" then
               return ....

           elif col="G" then
               return ....

           elif col="R" then
               return ....
          end if;
          end proc:

More controls could be necessary for a safer code (what if the denominators in the green and red cases are null?)

Download arctan.mw

Unfortunately I gave only Maple 2015 and WeibullPlot is not avaliable.

Look to this file and tell me if the plot I draw corresponds or not to what you would like to get?
I'm afraid I did not corerctly understood what is the issue you refer to. 
Here axis 1 is given 4 subticks and axis 2  3 subricks, you can see that the plot is correct.

Download WeibullPlot_2015.mw

but an interval.

Just do this to see what happens for a given y as x goes to 0.

limit(f(x, y), x=0)
                            /1\             /1\
                   -1 + sin| - | .. 1 + sin| - |
                            \y/             \y/

Maybe doing this will help you understand what happens

limit(f(1/x, 1/y), x=+infinity)
                   -1 + sin(y) .. 1 + sin(y)
limit(sin(x), x=+infinity)
                            -1 .. 1

You can also compute a directional limit 

limit(f(x, lambda*x), x=0)
                            -2 .. 2

In the attached file you will find an exploration of f(x,y) as you get closer and closer to the point (0, 0) (the value h of the slider means the neighborhood of (0, 0) has size 10^-h). 
limit.mw

A quick workaround is

plots:-implicitplot(g=1e-8, x = 0 .. 80, y = 0 .. 60, scaling = constrained, grid=[400, 400])

Nevertheless do not gorget that g=0 on the domain  33 <= x <= 58,  0, <= y <= 49.99226668
gnull.mw

(works only is the pattern to search is numeric, if not adaptations have to be done)

Download LookUp.mw

Note that I did not pay attention to performance and that the LookUp function may be inefficient for large matrices.

In case you would like to remove the columns/rows which contain the pattern to search, here is variant of the previous file LookUp_variant.mw

... here is a solution  I wonder what you mean by "The result will be a list of 56 models with 5 monomials."?)

It is efficient while the lengths of the two lists are not much larger:

restart:
with(combinat):
L1 := [x^2*y*alpha[1, 11], x*z^2*alpha[2, 15], y^2*z*alpha[3, 17] + y*z*alpha[3, 8]]:
         
L2 := [alpha[i, 0], alpha[i, 1]*x, alpha[i, 2]*y, alpha[i, 3]*z, alpha[i, 4]*x^2, alpha[i, 5]*y*x, alpha[i, 6]*z*x, alpha[i, 7]*y^2, alpha[i, 8]*z*y, alpha[i, 9]*z^2, alpha[i, 10]*x^3, alpha[i, 11]*y*x^2, alpha[i, 12]*z*x^2, alpha[i, 13]*y^2*x, alpha[i, 14]*z*y*x, alpha[i, 15]*z^2*x, alpha[i, 16]*y^3, alpha[i, 17]*z*y^2, alpha[i, 18]*z^2*y, alpha[i, 19]*z^3]:

#-----------------------------------------------------------------

# Build the cartesian product of the two lists

T:=cartprod([L1, L2]):
L1L2 := NULL:
while not T[finished] do 
  L1L2 := L1L2, T[nextvalue]() 
end do:
L1L2 := [L1L2]:

# Keep only the elements which do not have repeated alphas with the same second index.

map(u -> map2(op, 2, [select(type, indets(u), `indexed`)[]]), L1L2):
Keep := zip((u, v) -> if `not`(member(v[-1], {v[1..-2][]})) then add(u) end if, L1L2, %):

print~ (Keep): # to get Keep displayed one element per line

numelems(Keep)
                               56

A_solution.mw

If your matrix has the structure you give, that is +/- a given quantity at some positions and 0 elsewhere, here is a way=

Download For_instance.mw

For more complex matrices there are some ambiguities about what is the "common factor" to put outside of the matrix: For_instance_2.mw

restart

# list of small greek letters
seq(cat(`#mo("&#`, i, `;")`), i=945..969)



# A procedure to generate a random word of L letters
RandomGreekWord := proc(L)
  local r := rand(945..969);
  cat(
       `#mrow(`, 
       seq(cat(`mo("&#`, r(), `;"),`), i=1..L-1),
       cat(`mo("&#`, r(), `;"))`)
  )
end proc:

# an example
RandomGreekWord(10)

Download Greek.mw

Capital greek letters range from 913 to 937.
So all greek letters can be displayed by replacing

i=945..969

by 

i in GreekIndices

where 

GreekIndices := [$913..937, $945..969]