Using vpasolve to solve equations returns empty sym 0-by-1 although there is a solution

syms f a b
f = a*x^(-b);
%a=38.76;
%b=3.969;
%c=int(f,x,[3 Inf]);
%d=int(f,x,[4.5 Inf]);
%e=double(c)
%f=double(d)
eqns=[int(f,x,[3 Inf]) == 0.5003,int(f,x,[9/2 Inf])==0.1501];

vars = [a b];
[sola, solb] = vpasolve(eqns,vars)

I just don't understand why I get empty sym 0-by-1 as an answer. This upper fat tail of a Pareto distribution is supposed to meet the conditions of probability 0.5 and 0.15 for a~=38.76 and b~=3.969 integrating from 3 to infinity respectively from 4.5 to infinity. I'd appreciate very much any helpful suggestions!

Cheers, Sonja

%syms f a b
f = a*x^(-b);
a=38.76;
b=3.969;
c=int(f,x,[3 Inf]);
d=int(f,x,[4.5 Inf]);
e=double(c)
f=double(d)
%eqns=[int(f,x,[3 Inf]) == 1/2,int(f,x,[9/2 Inf])==15/100];

%vars = [b a];
%[solb, sola] = vpasolve(eqns,vars,[3.9,40])

Answers: 
e =

    0.5003


f =

    0.1501

Answer 1

[Sonja Glas]

I found the solution myself. The init_guess in vpasolve did somehow not help (maybe the expression was faulty?)

BUT assume and assumeAlso brought the solution:

syms f a b
f = a*x^(-b);
assume(b>3);
assumeAlso(a>35);
%a=38.76;
%b=3.969;
%F=int(f,x)
eqns=[int(f,x,[3 Inf]) == 1/2,int(f,x,[9/2 Inf])==15/100];

vars = [a b];
[sola, solb] = vpasolve(eqns,vars)

sola =

38.759579642286119630010875007585

solb =

3.9693622959161180276359691433877